Alban Leroyer
email: alban.leroyer@ec-nantes.fr

Michel Visonneau

Sophie Barré

J.-M Kobus

Jean-Michel Kobus

Experimental and numerical investigations of the fl ow around an oar blade

Keywords: Rowing, Oar blade, Free-surface fl ow, Navier-Stokes solver

framework is presented and the Navier-Stokes solver and methods for handling multifl uid fl ows and moving bodies are described. Lastly, numerical results are compared with experimental data, highlighting an encouraging agreement and proving the relevance and the complementarity of both approaches.

Introduction

Fluid mechanics is a scientifi c fi eld that plays a key role in nautical sports, mostly in complex situations. This is particularly true for rowing, for which fl uid fl ow knowledge is required twice: once around the boat hull and again around the oar blades. Moreover, all the phenomena involved in rowing (the movement of rowers, the forces acting on the oars, and the behaviour and resistance of the hull) strongly interact. Consequently, the only means of analysing and understanding rowing is to use dynamic simulations. Therefore, since 1998 the Fluid Mechanics Laboratory of Ecole Centrale de Nantes (LMF) has developed a simulator for a complete boatoarsmen-oars system. This tool constitutes a means of forging a permanent and operational synthesis of the research and studies on rowing. The accuracy of simulations evolves according to the progress in results and validations. The relevance and the applicability of our simulator required a thorough study of the phenomena involved in oar propulsion. Large variations in speed added to the heave and pitch of the boat induced by the movement of oarsmen during the rowing stroke, make the fl ow around rowing boats diffi cult to investigate. Abstract This article aims at verifying the capabilities of a Reynolds-Averaged Navier-Stokes Equations (RANSE) solver (ISIS-CFD, developed at the Fluid Mechanics Laboratory of Ecole Centrale de Nantes [LMF]) to accurately compute the fl ow around an oar blade and to deduce the forces on it and other quantities such as effi ciency. This solver is structurally capable of computing the fl ow around any blade shape for any movement in six degrees of freedom, both when the blade pierces the free surface of the water and when it does not. To attempt a fi rst validation, a computation was performed for a simplifi ed case chosen among those for which experimental results are available at LMF. If results prove satisfactory for a simplifi ed blade shape and for a movement that respects the main characteristics of blade kinematics, then the solver could be used for real oars and more realistic kinematics. First, the experimental setup is considered, and the objectives, methodologies, and procedures are elucidated. The choice of the test case for numerical validation is explained, i.e., a plane rectangular blade with a constant immersion and a specifi ed movement deduced from analogy with tests on propellers. Next, the numerical However, the fl ow around an oar blade appears far more complex and has been far less studied.

Most authors adopt the very simple model suggested by Wellicome [START_REF] Wellicome | Some hydrodynamic aspects of rowing[END_REF] to evaluate the force on oar blades and to develop considerations on optimisation of oar setting and rowing style. This model supposes that the force is perpendicular to the blade and proportional to the square of the normal component of the instantaneous absolute velocity calculated at the centre of the blade. This quasi-static model is dimensionally consistent, but the coeffi cient of proportionality is never precisely specifi ed. Some computational fl uid dynamics (CFD) studies have been carried out using Videv [START_REF] Videv | Unsteady viscous fl ow simulation around the blade of a rowing boat[END_REF] to calculate this coeffi cient in a two-dimensional case without a free surface, but this simplifi ed confi guration is too far from the specifi city of fl ow around oar blades to be helpful.

To improve the evaluation of efforts on blades induced by such a complex fl ow, only an experimental approach seemed to be feasible until a few years ago. To perform such experiments, the LMF built a sophisticated device for testing blades at reduced scale and also for testing real scull oars in its towing tank. Nevertheless, during the past decade, thanks to the growth of storage capacity and computing power, CFD tools have become no longer limited to simple physical problems and have signifi cantly enlarged their fi eld of applications. The in-house Reynolds-Averaged Navier-Stokes Equations (RANSE) solver ISIS-CFD integrates new physical features to deal with increasingly realistic applications, such as moving bodies and complex free-surface fl ows. Thus, it seemed to us that the time had come to compare the results of CFD simulations with experimental results on such a complex confi guration. As a result, initial simulations have been recently performed. To explain our approach to the problem and to comment on the results are thus the main goals of the present article, which is devoted to the experimental and numerical aspects of oar blade hydrodynamics.

First, Sect. 2 describes the experimental approach, the procedure, and the aims of such a study. The experimental confi guration chosen for comparison with the numerical results is one of those tested during a testing campaign, the aim of which was to build a model of hydrodynamic forces on a blade as a function of the movement in the water and the geometric characteristics. To discuss this modelling is not the subject of this article; nor is it intended to make any conclusions about real rowing styles. Here, the test bench and the associated procedures are described, so that the reader can evaluate the degree of confi dence which can be afforded to the experimental results. The blade chosen for the study is the reference blade of systematic tests, that is to say, a plane rectangular blade with a rigid shaft moving with analytical kinematics defi ned by only two parameters. The reasons for this selected confi guration for numerical validation conclude Sect. 2.

In Sect. 3, the numerical methods are considered and the main features of the RANSE solver are given. Special attention is paid to the numerical treatment of the moving grid in the Navier-Stokes equations. This point needs to be addressed as soon as body motion is encountered. Next, the simulations performed for this initial study on an oar blade are detailed in terms of meshes, boundary conditions, law of motion. The time evolution of the free surface is lastly described. In Sect. 4, comparisons between available experimental measurements and simulations are shown. Forces and torque are considered fi rst. A comparison of experimental and numerical evaluation of effi ciency completes the discussion. Lastly, conclusions highlight the complementarity of both approaches and the perspectives of such research.

Experimental approach

General methodology of the test campaigns

Instrumented boats are now commonly used to measure the motions and forces acting on oar blades. Such tests have already been performed, in particular at LMF. They have the advantage of working directly on a realistic confi guration, but the complexity of in-situ measurements, which involve numerous and sometimes not precisely controlled parameters, does not allow adequate precision and repeatability to build or validate models of forces on blades. For example, it is technically difficult to measure the hydrodynamic force on a blade, and the position (or velocity) of blades is not precisely known because of oar fl exibility and variable immersion. Thus, in-situ measurements have been preferentially used until now for comparative analysis and for global validation of rowing simulations. Another approach is to reproduce as well as possible the movements of blades in a towing tank. This technique permits the separation of the propulsive device (the oar) from the motor (the oarsman) and from the boat, as it is done for propellers. Such a system was designed at LMF 3 for testing real oars or models of blades at a reduced scale with a rigid shaft. The system generates a simplifi ed rowing stroke but with good control of the parameters and with accurate reproducibility. In order to avoid a complex mechanism, the blades always remain in the water. To limit the consequences of this feature, it is possible to impose a movement (denoted the neutral movement) calculated to minimise the drag on the blade before the catch angle of the real stroke. This technique introduces supplementary parameters such as the neutral motion duration and the transition phase duration (the catch). In order to perform systematic tests, we also defi ned simplifi ed movements with only two kinematic parameters specifi ed under efficiency considerations, as described later.

Experimental device description

The device is a six-component dynamometer equipped with a servomotor and a mechanism which reproduces accurately and repeatably the rotation of the oar blade in the water (see Fig. 1). The servomotor is a brushless motor with high specifi c power. Reduction of the rotation is done with a gear box associated with two pulleys and a belt. The motor has a numerical control system which ensures that the motor conforms to a given displacement (or speed) reference signal. The dynamometer allows the servomotor and the transmission mechanism to be fi xed directly on to the measuring device. This arrangement eliminates the infl uence of interior transmission efforts and allows the measurement of the hydrodynamics force on the blade, corrupted only by some parasitic effects. The device was fi xed to the towing tank carriage and the speed of translation was constant during the stroke. Immersion adjustment was set by shifting the dynamometer support up and down. The vibrations of the carriage excite the natural frequencies of the dynamometer. From analysis of some in-situ measurements on a skiff, we deduced that the amplitude spectrum of the forces on oars has no signifi cant value higher than 7 Hz, so we designed the device to have natural frequencies higher than 15 Hz. To achieve these favourable conditions, the measuring part of the dynamometer and the rotary arms were built from carbon fi bre. The test bench measured forces and moments acting on the oar in conjunction with the characteristics of the predefi ned movements, but deformation of the oar shaft and defl ection of the blade were not measured. As the position of the blade in the water has to be precisely known for numerical validation, we chose to use results obtained with the bench equipped with a rotating arm suitable for blades of reduced scale (see Figs. 1 and2). Along with the arm, the blade shaft is very short and was made of monolithic carbon to obtain great rigidity. It was fastened on the rotary arm using a device that makes it possible to adjust the external lever, the inclination of the oar shaft and the angle of the blade in the vertical direction.

Experimental procedures

Before each run, the centre of gravity was adjusted to lie on the rotation axis and the inertial characteristics were measured. The inertial effects consisted of moments only, because the translation speed was constant. We also measured the aerodynamic forces on the device. The numerical fi le of reference signals used by the control system was calculated with a specifi c program that computed the different phases of the blade rotation during the experiment. A run proceeded in four steps. The oar blade was initially oriented at 0°. The carriage was then accelerated to the desired speed, and then ran at this constant speed. Then, the blade followed the specifi ed trajectory up to an angle of 180° before the carriage stopped. The measurement processing was done as follows. The device produced eight analogue voltage outputs: six outputs came from force transducers in the dynamometer and the two other outputs were the angular position and velocity of the servomotor. These signals were sampled and stored. The six components of forces and moments in the dynamometer frame of reference were calculated at each sample time using a calibration and coupling matrix. Next they were fi ltered to eliminate noise and forces induced by vibrations in the dynamometer. The angular position and velocity were also sampled and fi ltered. Finally, we introduced corrections to eliminate the inertial effects induced by experimental proce- Sampling and measurement acquisition. For measurements acquisition, we use an A/D converter card controlled by an advanced program. This program permits the defi nition of acquisition sequences and stored them. The sampling frequency was chosen to be 500 Hz to avoid using an anti-aliasing fi lter.

Filtering. The fi ltering method consisted of a lowpass windowing spectrum followed by a reconstitution of the signal by applying an inverse Fourier transform. This technique did not induce any phase shift between measured signals. A band of frequency between 7 and 15 Hz having very low spectral power allowed low-pass fi ltering to be performed with good precision and with few induced perturbations in the pass band.

Pre-processing. This procedure consisted fi rst in ordering and storing raw data and fi ltered data, and second in correcting data using inertial and aerodynamic effects.

Measurement exploitation. First, the physical parameters used for analysis and modelling were computed, and second, kinematic data and forces measured in the dynamometer axis were used to calculate the forces in a different axis (e.g., the blade axis or the relative incidence of water). Finally, hydrodynamic coeffi cients, the instantaneous driving force, the effective power, and the instantaneous and global effi ciencies were evaluated.

Blade kinematics and effi ciency considerations

To begin with, some notations were defi ned as follows and are represented in Figs. 234.

Geometry and kinematics

• (x c , y c , z c ): frame of reference linked to the carriage (or boat), with z upward and the motion of the carriage directed towards -x c

• ( t , n, z): frame of reference linked to the blade, with t the horizontal vector tangent to the blade directed towards the rotation axis • L e : outboard lever (>0), i.e., the distance between the rotation axis (O, z) and I, the centre of the blade surface

• q = (xt c , ): the angle between the shaft and the translational boat velocity • b: the angle between hydrodynamic force F h and n, the vector perpendicular to the blade

• V b = | V b | = -V b • x c : the boat velocity (>0).

Components of horizontal hydrodynamic force F h and the torque in different axes of reference (see Fig. 4)

• F x , F y : in the dynamometer (or boat) axis • F t , F n : in the shaft axis, • M z : the component of the hydrodynamic torque applied about (O, z).

The instantaneous effi ciency h is defi ned as the ratio of the useful power to the delivered power:

η θ = FV M xb z ɺ (1)
Let us denote the centre of the blade I, and the point where the resultant fl uid force F h acts at each instant P (see Fig. 4). Then, Eq. 1 can be formulated as:

η θβ βθ θβ = -⋅ + () ⋅ + ( ) ⋅⋅() = + () ⋅ - FV OI d F V L hb ip h b e ɺ sin cos sin d d ip ( )⋅ () cos , βθ ɺ (2)
where d ip is the algebraic distance between I and P projected on the oriented shaft axis (O, t ). Note that M z = -(L ed ip ).F n = (L ed ip ).(F x sin q -F y cos q). So, we can deduce:

By extending Eq. 2, we can thus obtain:

η θ θ β θ ηη η =+     -       =⋅⋅ V L dL b e ip e fp sin tan tan ɺ 1 1 1 0 (4)
This defi nes three components of effi ciency:

ηθ θ η β θ η 0 1 11 = ( ) =+ ( ) =-( ) VL dL be f pi p e sin , tant a n, ɺ (5)
The effi ciency component h p , related to the centre of pressure, is close to unity since d ip cannot be greater than half the length of the blade, i.e., quite small relative to L e . The component h f depends on the orientation b of the hydrodynamic force F h relative to the blade perpendicular n. It is always unity when q = p/2. During the angular range of a real stroke (where q can vary from about 40° through to 130°), b remains small and h f is close to unity. Thus, the component h 0 appears to be the most signifi cant component of the effi ciency. Its value is close to that of h when F x and M z are maximum and it depends only on kinematics and on the characteristic length of the oar L e .

By analogy with propellers, note that h 0 can be regarded as a kind of advance number that is representative of the incidence history. More precisely, if V I = V n n + V t t is the velocity of the centre of the blade, we can defi ne the pseudo-incidence of water, denoted i, as the angle of the opposite of V I relative to the tangent vector t . Then we can write:

tan sin cos tan sin tan i LV V L V eb b e b = - =-       =-    ɺɺ θθ θ θ θ θ θ η 1 1 1 0     (6)
The history of the pseudo-incidence i depends only on the angular position of blade and on h 0 . It can be noted here that during the possible phase of neutral movement (which corresponds to h 0 = 1 for a q value less than the catch angle), the pseudo-incidence remains equal to zero.

To defi ne blade movements for systematic tests, we cannot reproduce real measured oar movements because these measurements are performed with oarsmen having different rowing styles and are not linked with reliable data on instantaneous boat speed. Furthermore, the high number of parameters that defi ne this kind of movement makes the exploitation of test results diffi cult. We observed great disparities in the forces and moments curves according to chosen movement parameters. Con-sequently, to perform systematic tests, we defi ned movements with two main parameters according to kinematics and energetic considerations. Wellicome, 1 and Pope [START_REF] Pope | On the dynamics of men and boat and oars[END_REF] after him, showed that the global effi ciency of a stroke has a maximum value when the instantaneous effi ciency is constant during the stroke. Wellicome emphasized that Betz used this argument originally in determining the optimum distribution of pitch for a propeller. So we decided to test oar blades with movements that make the component h 0 constant. This is the simplest way of keeping effi ciency h almost constant during a wide part of the stroke (including the angular range of a real rowing stroke), since h 0 depends only on kinematics and geometry. Modelling of efforts will be done in this condition. This choice seems to be restrictive, but elite rowers search for optimum effi ciency for a given power and their rowing style is probably close to one of constant effi ciency. Then, if h 0 remains constant during the stroke, the defi nition of h 0 in Eq. 5 imposes the following condition:

∀= = θθ θ θ ,s i n , ɺɺ KK with max (7)
Thus, we have:

η 0 = V LK b e (8)
Consequently, systematic tests were performed for different pairs of parameters (K, h 0 ), the K value varying from 3 to 4 rad•s -1 and h 0 varying from 0.59 to 0.8. The boat (carriage) velocity was deduced from Eq. 8. By integrating Eq. 7, we obtain:

θ a tK t t ()=-( ) ( ) [ ] 2 0 atan exp , ( 9
)
with t 0 the time where q a = p/2. For t = 0, q a (0) ≈ 2. exp(-Kt 0 ) is very small but not zero (typically, q a (0) = 0.01°, q . a (0) ≈ 7×10 -4 rad•s -1 ). Between q = 0 and q a (0), a brief linear junction with a slope value of K sin(q a (0)) is applied to reduce the velocity discontinuity at q a (0). The time offset ∆t = 1/K introduced by the linear junction is taken into account to generate the reference signal of rotation. When h 0 is constant (h 0 = Cst), according to Eq. 6 and Eq. 9, we see that the pseudo-incidence history i(t) is completely defi ned by the characteristic effi ciency (or advance number) h 0 and K = q . max . Figure 5 gives an example of the temporal evolution of q for K = 3 (t 0 = 3.1084 s). Figure 6 shows the deviation between the analytical angle q a and the measured value q m . Because of the a posteriori recalibration of the analytical law q a with measurements, which defi nes a corrected time t 0 and ensures that q m (t 0 ) = q a (t 0 ) = p/2, the error q mq a is then approximately zero around time t 0 . As we can see, the maximum deviation is quite small (about 0.2°). Analysis of the slope shows that the maximum error in terms of q .

does not exceed 0.7°s -1 in a range where q . is always greater than 70°s -1 , giving a relative deviation of less than 1%. Thus, the analytical law has been used to prescribe the motion in CFD simulations instead of the measured data.

Choice of blade and experimental confi guration for CFD comparison

Systematic tests were performed with six reduced-scale blades. In order to simplify the meshing of the blade, we chose the reference blade, that is to say, a planar rectangular blade (see Fig. 2), to compare measurements with numerical results. This choice does not alter the nature of the fl ow or the generality of comparisons. The parameters of the chosen movement were h 0 = 0.5887 and K = q . max = 3 rad•s -1 . According to Eq. 8, the velocity of car-riage was 2.31 ms -1 . The time evolution of q is shown in Fig. 5 and the blade trajectory in Fig. 7. A particularly low effi ciency with the blade close to the free surface was chosen because the fl ow is more perturbed than with highly effi cient motion. With the aim being to validate the numerical approach, it seemed to be interesting to evaluate the capabilities of ISIS-CFD in calculating a fl ow with such severe conditions. If the agreement with experimental data turned out to be good for this case, we could be hopeful that the code, validated in such a complex confi guration, would provide reasonable results when simpler confi gurations characterized with more realistic higher effi ciency are studied.

Summary of hypothesis and simplifi cations

-The oar blade is a plane rectangle with a thickness of 3.5 mm (see Fig. 2). -The area of the blade corresponds to a scull oar at scale 0.7 with the same aspect ratio. -The blade chord is aligned with the shaft.

-The blade is vertical and its upper part remains at the level of the still free surface. -The deformation of the shaft and blade is supposed to be negligible. -The translation has a constant velocity.

-The movement is specifi ed with a constant characteristic effi ciency h 0 (or advance number), without a neutral phase at the beginning of the stroke.

It can be noted here that the specifi city of the fl ow is mainly due to the unsteady rotation of the oar combined with the translation movement of the boat. It gives typical trajectories as shown in Fig. 7. As a consequence, despite the simplifi cations explained above, the fl ow obtained around the blade is representative of the greater part of the physical phenomena that occur during a real oar stroke. This is a major point, since the objective is to validate numerical modelling of such a fl ow.

3 Numerical approach

Flow solver

The ISIS-CFD solver, developed by Equipe de Modélisation Numérique, the CFD Department of LMF, uses the incompressible unsteady Reynolds-averaged Navier-Stokes equations (URANS). The solver uses the fi nitevolume method to build the spatial discretization of the transport equations. The face-based method is generalized to unstructured meshes for which non-overlapping control volumes are bounded by an arbitrary number of constitutive faces. The velocity fi eld is obtained from the momentum conservation equations and the pressure fi eld is extracted from the mass conservation constraint, or continuity equation, transformed into a pressure equation. In the case of turbulent fl ows, additional transport equations for the modelled variables are solved in a form similar to the momentum equations and they can be discretized and solved using the same principles. Several turbulence models, ranging from the oneequation Spalart-Allmaras model, [START_REF] Spalart | A one-equation turbulence model for aerodynamic fl ows[END_REF] two-equation kw closures, [START_REF] Menter | Zonal two-equation k -w turbulence models for aerodynamic fl ows[END_REF] to a full-stress transport R ijw model, [START_REF] Duvigneau | On the role played by turbulence closures for hull shape optimization at model and full scale[END_REF] are implemented in the fl ow solver to take into account turbulence phenomena. Incompressible and non-miscible fl ow phases are modelled through the use of conservation equations for each volume fraction (or concentration) c i of each phase i. Considering the incompressible fl ow of a viscous fl uid under isothermal conditions, mass, momentum, and volume fraction conservation equations can be written as follow (using the generalized form of Gauss' theorem):

U ndS ⋅= ∫ 0, S (10a) δ δ ρρ ρ t UdV U UU ndS gP d V d t VS VS ∫∫ ∫∫ +-( ) ⋅ =-∇ ( ) ++ ( ) ⋅ ޔޔ n ndS, (10b) δ δt cdV c UU ndS ii d VS ∫∫ +-( ) ⋅= 0, ( 10c
)
where V is the domain of interest, or control volume, bounded by the closed surface S moving at velocity U d with a unit normal vector n directed outward. U and P represent, respectively, the velocity and the pressure. ޔ and ޔ t refer to the viscous and Reynolds stress tensors, whereas g is the gravity vector. The time derivative fol-lowing the moving grid is written d/dt. Whereas ޔ t is determined according to the turbulence model used, ޔ follows the classical relation of a Newtonian fl uid for incompressible fl ows. The effective fl ow physical properties (viscosity m and density r) are obtained from each of the phase properties (m i and r i ) using the following constitutive relations:

ρρ µµ == = ∑∑ ∑ cc c ii i ii i i i ,, 1 (11)
Except for the convection terms and volumetric mass fl uxes, the interfacial quantities q f are rebuilt linearly from the quantities themselves and their available cellcentred gradients. Special attention has to be paid to face reconstructions of the volume fraction c i . The challenge raised by the discretization of a transport equation for the concentration is the accurate modelling of the contact discontinuity, i.e., the free surface. In order to assume face-bounded reconstructions and to avoid unrealistic oscillations, the search for an acceptable compromise between accuracy and boundedness of c i ∈ [0, 1] is a key point. 9,10 Moreover, the method must also preserve the sharpness of the interface through the transport equation (Eq. 10c). These requirements are fulfi lled by the inter-gamma differencing scheme (IGDS), [START_REF] Jasak | Interface tracking capabilities of the inter-gamma differencing scheme[END_REF] which introduces downwind differencing since compressive characteristics are required for accurate interface capturing. Through a normalized variable diagram (NVD) analysis, 12 this scheme enforces local monotonicity and the convection boundedness criterion (CBC). [START_REF] Gaskell | Curvature compensated convective transport: Smart, a new boundedness-preserving transport algorithm[END_REF] The main disadvantage of the IGDS scheme (and also of other compressive schemes) is a Courant number limitation, Co < 0.3 in multidimensional cases, known as the Courant-Friedrich-Levy (CFL) condition. The Courant number of any face is defi ned as follows:

Co t V =∆ F ( 12
)
where F is the total (positive) velocity fl ux through the considered face, V is the volume of the upwind cell, and ∆t is the global time step of the temporal discretization.

Therefore, the discretized time step has to be small enough to fulfi l the CFL condition. It has been shown that the role played by the compressive property of the IGDS is fundamental to achieve a reliable simulation of the free surface. [START_REF] Queutey | An interface-capturing method for free-surface hydrodynamic fl ows[END_REF]

Body motion

To implement a body's motion in a fl ow solver, the mesh must be adapted to the new position of the body over time. In order to keep an appropriate grid, three complementary methods have been integrated (for further details, see Leroyer 15 ):

-the pseudo-structure regridding procedure, -the rigid transformation of the mesh, -the analytical weighted regridding approach.

When the grid is moving, the so-called space conservation law must be satisfi ed:

δ δt dV U ndS d VS ∫∫ -⋅ = 0 (13)
The mesh mobility is then taken into account in the equations of conservation with the grid velocity U d , or more precisely with the grid displacement velocity fl ux through each face S, denoted by F S U d . It is defi ned by:

F S S Ud d U ndS =⋅ ∫ (14)
For deformation techniques (pseudo-structure analogy and weighted regridding), this quantity is obtained by computing the exact volumes swept by cell faces during the time steps, which ensures that the discrete space conservation law is exactly satisfi ed. [START_REF] Leroyer | Fluid/motion interaction for solid and fl exible bodies by resolution of the Navier-Stokes equations[END_REF][START_REF] Ferziger | Computational methods for fl uid dynamics[END_REF] If a rigid transformation is also applied, its contribution to the displacement velocity fl uxes can be exactly computed by using the parameters of the kinematic screw {Ω, U } O 1 of the rigid transformation. Indeed, for each point M geometrically linked to this mesh, the velocity induced by the rigid displacement is defi ned by:

UM UO O M dd () = () +⊗ 11 Ω
This property is used to generate a direct and exact calculation of the grid displacement velocity fl ux through each face S using the following:

F F S SS Ud d f r d U ndS U F S FM ndS =⋅ = () ⋅+ ⋅ ⊗ ∫∫ Ω ' (15)
F denotes the centre of the face and S f the face vector. For a planar face (composed of three nodes), the defi nition of F is such that the integral of the term F r is strictly equal to zero. As far as non-planar faces are concerned (composed of more than three nodes), they are split into triangular faces, in which F is their common summit (see Fig. 8). By summing the contribution of each triangle making up the face, the term F r can be evaluated. In this way, the calculation of F S U d is exact, whatever the face considered. As the motion is a rigid displacement, we have for any point M geometrically linked to the mesh:

div U d ( ) = 0 (16)
Applying the divergence theorem after integrating Eq. 16 yields: [START_REF] Leroyer | Numerical methods for RANSE simulations of a self-propelled fi shlike body[END_REF] For simulations dealing with an isolated solid body in an infi nite domain (as is the case here), only a rigid transformation is employed. This technique has no major effect on the computing time and enables motion of unlimited amplitudes in the physical space. Since there is no deformation, each cell conserves its volume. Therefore, the unsteady term dV/dt is identically null for rigid displacement. Considering Eq. 17, the space conservation law is then exactly satisfi ed.

div UMd V UMndS dd U d () ( ) = = () ⋅= ∫ ∫ 0 ∂V V F S S faces ∑ = 0
In the case of prescribed motions, at every time step, bodies are fi rst displaced, then the mesh is regridded, and lastly the fl ow is solved. Hence, there is no real coupling, and this is because of the lack of fl uid feedback on the body position. Since the motion is imposed here, the question of the resolution of the motion and of the coupling with the fl ow is not considered. For solved motions, further information is given by Leroyer [START_REF] Leroyer | Fluid/motion interaction for solid and fl exible bodies by resolution of the Navier-Stokes equations[END_REF] and Leroyer and Visonneau. [START_REF] Leroyer | Numerical methods for RANSE simulations of a self-propelled fi shlike body[END_REF]

Characteristics of simulations

Three simulations were performed; their characteristics are summarized in Table 1. The fi rst simulation was based on a coarse mesh of 120 000 cells with an adaptive time step law to ensure the compressive property of the

S f

Fig. 8. Splitting of a non-planar face volume fraction scheme, whatever the fl ow velocities. [START_REF] Hay | H-adaptive Navier-Stokes simulations of free-surface fl ows around moving bodies[END_REF] The second was performed on a much fi ner grid (M14e5) of 1 400 000 cells, but with a 0.005 s constant time step, which lead to Courant numbers locally larger than 0.3. The third simulation also used this fi ne grid, but with an adaptive time law. An insight into the mesh fi neness can be seen in Fig. 9, which represents the surface grid on the blade. It can be noted on Fig. 9b that even if it contains 1 400 000 cells, the fi nest grid M14e5 has not so many small cells because the whole domain around the blade has to conserve small cells to capture the free surface precisely. This is because the free surface moves along the whole surface of the blade during the motion. Therefore, in quite a large area around the blade where the free surface is largely in motion, cells of small size have to be kept. In fact, without an automatic local adaptive mesh refi nement, as has been detailed by Hay et al., [START_REF] Hay | H-adaptive Navier-Stokes simulations of free-surface fl ows around moving bodies[END_REF] this feature is necessary to maintain a well-defi ned (i.e., not too diffused) interface during the whole stroke. A global view of the boundaries of the fl uid domain (except the upper boundary) is shown on Fig. 10. The surface of the blade can be observed in the middle of the fl uid domain. In fact, the blade is not positioned exactly in the middle because the motion is not symmetric. Figure 11a illustrates this point: it shows a horizontal slice of the mesh (fl uid domain is vertically extruded). Figure 11b is a close-up of Fig. 11a around the blade. The boundary conditions are imposed as follows: on the top, a hydrostatic pressure is imposed, whereas a slip condition is required on the bottom of the domain (a hydrostatic pressure could also have been imposed here).

In the lateral regions, the velocity is imposed with its far-fi eld value. For these fi rst simulations, only a slip condition was tested on the blade, neglecting the turbulence effects around it. In this study, the Earth was chosen as the Galilean frame of reference. Hence, the translation of the carriage and the rotation of the blade are prescribed as an imposed motion of the blade (see Fig. 7) and the far-fi eld velocity is then supposed to be at rest. Another possibility would have been to work in a Galilean frame linked to the carriage and to impose on the blade only the motion of rotation (and also a velocity equal to V x b for the far fi eld). But in the fi rst case, the acceleration phase of the carriage can be reproduced more easily. Here, before launching the motion of rotation, a transient state for the motion of translation was used to raise smoothly the V x component from zero up to the constant speed of the carriage. As detailed in Sect. 3.2, a rigid transformation of the mesh was favourably used: cells follow exactly the motion of the blade without deforming; the mesh quality is then conserved in time. The simulations exhibited the very complex evolution of the free surface. Figures 12-17 represent the free surface at different times in the blade frame. The (a) parts of the fi gures are viewed from the side away from the shaft, looking towards t and blanking values on the t axis lower than -1.458 m (the origin being set at point O, the centre of rotation of the oar, see Fig. 3). The centre of view of the (b) parts of the fi gures is located on the shaft side, looking towards the blade, and blanking values greater than -1.158 m, i.e., a quarter of the blade spanwise. At the beginning of the rotation, it can be observed in Fig. 12 that there is an elevation of the free surface just behind the infl ow side of the blade, and consequently a fl ow passing over the upper part of the blade. Simultaneously, the low-pressure level generated on the other side of the blade leads to a bowl-shaped cavity that increases in size. The cavity and the elevation become more intensifi ed as time (and so q) increases and the fl ow above the blade gradually takes the shape of a conical breaking wave (see Figs. 13, 14). It is noted that this conical wave points towards the external tip of the blade. The deformation is thus more accentuated 15), generating a huge ventilated cavity. We can see here that the free surface seems to be more disturbed compared with a real rowing stoke. Videos of experiments (with a camera fi xed to the shaft) tend to confi rm this observation. This phenomenon is not so surprising since the effi ciency h 0 used here (=0.5887) is lower than that of real strokes, thus generating a more violent fl ow. In fact, the higher the effi ciency, the less perturbed the fl ow is. The cavity is then convected (see Fig. 16) while the tip of the conical-shaped wave reverses its position with respect to the blade (now pointing towards the shaft). At the end of the motion (see Fig. First, looking at the period before the rotation is significant (and at the end of the stroke), we can see that the simulation gives lower drag in the x-direction. This is not surprising, since the motion of the blade is like that of a fl at plane moving straight ahead (due to the translation of the carriage) in this part of the stroke. At this stage, the viscous drag acting on the blade is much larger than the pressure drag. But, because a slip condition is imposed on the blade, the viscous drag is not taken into account (calculated forces are mainly due to the pressure contribution), which explains the observed discrepancy.

Focusing on the F x component, which represents the propulsive force, it can be seen that the Sim3 simulation gives better trends compared to the two others. This is very encouraging because Sim3 is the simulation with the fi nest grid using an adaptive time step leading to activation of the compressive properties for the numerical scheme applied to the discretization of the volume fraction transport equation. Even if grid independence is unlikely to be achieved, the time evolution of forces is quite encouraging. In fact, the maximum value of F x is closer to the experimental value and is quite well reduced compared to Sim1 and Sim2. The time evolution is also better reproduced, even if it is not in total agreement. But Sim3 is the only simulation able to capture an infl exion at about t = 2.85 s. Experimentally, this phenomena seemed to be found earlier at t = 2.65 s. This time discrepancy might be explained by the turbulence effect of the dynamics of the fl ow, especially at those angles where a dynamic stall occurs. This latter point is perhaps even more visible in Fig. 19a, in which the longitudinal forces are plotted as a function of the angle: an experimental break is visible at about 25°, whereas the same phenomenon seems to appear later numerically (at about 45°). Similar conclusions can be drawn on the F y component (see Figs. 18b, 19b), for which Sim3 gives better results. However, there remain differences for higher values of F y obtained between t = 2.65 s and t = 2.85 s, i.e., between q = 25° and q = 50°. Surprisingly, it is noted that during the last part of the stroke, the evolu- tion is very well captured, especially when the forces stop decreasing at about t = 3.3 s. Comparison of the z-axis torque shows similar behaviour. Fairly good agreement is observed, apart from the angular zone between 25° and 50°. Here again, it seems that the break observed in the measured increase of M z at q = 25° and t = 2.65 s occurs about 0.2 s later in the simulation. In conclusion, we can say that the phenomenon observed experimentally at t = 2.65 s, (which results in a break on the evolution of forces F x and F y ), is not yet accurately captured by the simulations. In fact, it seems to occur later numerically (for the more precise simulation Sim3). These differences could perhaps be explained by too coarse a grid even on the M14e5 mesh and also by the infl uence of the turbulence effect on pressure and then on the free-surface elevation, especially when a dynamic stall occurs. It can also be noted in Fig. 6 that, after t = 2.6 s (up to t = 3.0 s), the difference between the analytical law (q a ) prescribed by the simulation and the real measured angle (q m ) varies rapidly. It locally leads to a difference of about 0.4 ° s -1 for the angular velocity, which is surely too small to explain the observed gap between the measured and simulated forces. due to h f , i.e., the orientation of the hydrodynamics with respect to the perpendicular direction of the blade. As Eq. 5 shows, this orientation has all the more infl uence on h as q deviates from p/2. For the range of angles corresponding to the fi rst part of the stroke (just after the catch), b is negative, i.e., the orientation of the force is unfavourable for propulsion. The propulsion is then lower than if the force were directed along n. The difference between the simulation and experiment here is found to result from quite a large discrepancy in h f in the fi rst part of the stroke, whereas the evolution of h p is in good agreement with measurements. As noted previously, the position of point P is actually well captured by the simulation. It is also seen that this coeffi cient h p remains close to 1 during a large part of the stroke (h p ∈ [0.94, 1.01] for q ∈ [40°, 140°]).

Figure 23 fi nally shows that h 0 can actually be considered the dominant parameter for such a stroke between q = 40° and q = 140°, corresponding to the effective propulsive phase for oarsmen. Here again, quite a large gap between the experimental and numerical efficiency h is visible in the fi rst part of the stroke (mainly due to the component h f , as remarked previously). For angles less than about 20°, the difference is quite large, because the viscous contribution to F x (which cannot be captured by the simulation since a slip condition is imposed on the blade) remains quite large compared to the pressure effect. In contrast, it does not infl uence F y similarly since the latter is mainly produced by the pressure effect when q is small. As a result, the angle b (and so h f ) cannot be precisely evaluated in this case.

Conclusions

The fl ow around an oar blade is very complex, especially because a human is part of the system. From the hydrodynamics point of view, it is an unsteady 3-D fl ow including violent free-surface motion and overturning. From the kinematics point of view, the motion of the blade results in a variable translation of the boat coupled to a rotation of the shaft and a vertical movement of the blade (even if in the present analysis, the latter is not taken into account). To investigate the fl ow and to improve the modelling of the forces acting on a blade (as part of a simulator for a complete boatoarsmen-oars system), a device reproducing a simplifi ed rowing stroke was designed and built. The chosen kinematics are not excessively restrictive and give a reliable database to validate the numerical approach. The numerical simulations presented here demonstrate the capacity to deal with such a complex fl ow, especially since the present work uses a low-effi ciency movement that produces a more violent fl ow than that of a real rowing stroke. Even if grid convergence is not obtained and turbulence is not yet taken into account, comparisons with experimental data are very encouraging and give confi dence in the capacity of the numerical method to simulate accurately such a fl ow pattern. The fact that the elevation and the shape of the free surface extensively change during the stroke makes the task more diffi cult. As a matter of fact, the whole part of the domain in which the interface takes place has to be gridded with fi ne cells. Without a local adaptive grid procedure, meshes with a huge number of cells are needed and simulations are hampered by a lot of useless fi ne cells. As has already been found by Hay et al., [START_REF] Hay | H-adaptive Navier-Stokes simulations of free-surface fl ows around moving bodies[END_REF][START_REF] Hay | Computation of free-surface fl ows with local mesh adaptation[END_REF] the development of a fully parallelized local mesh refi nement procedure with dynamic load balancing should decrease dramatically the processing time, while keeping the generality of the approach. With this technique, accurate simulations might be performed within a reasonable time. When this technique is operational and the processing time far shorter, new simulations will be performed for other experimental data of an oar with a rigid shaft. Comparisons with available data from real oars realized with the experimental device could also be added. Then, it will be possible to study the infl uence of oar deformation. Once simulations are validated with these towing tank tests, numerical simulations for a realistic confi guration will be relevant and will become far easier than experimental approaches. Indeed, the study of complex movements with variable immersion (with catch and fi nish) and with a variable speed of translation taken into account will become much easier than carrying out tests experimentally. The last step will then be validation by real rowing tests, which we have been carrying out for several years with instrumented oars and measurement equipment on the boat to record the oar and oarsman kinematics. With the numerical results of parametric studies on oars (e.g., stiffness, lever length, and shape of the blade) and movement characteristics, we will be able to replace empirical indications by objective and unbiased criteria not only for oars but also for the rowing style. Furthermore, a detailed analysis and comparison of the fl ow topology for some selected and typical results could then be started to precisely understand the physical mechanism involved in such a complex confi guration. The knowledge gathered in these fi elds and from such models as a result of these studies will be then progressively integrated into the global simulator.

Fig. 1 .Fig. 2 .

 12 Fig. 1. Photograph of the rowing device

Fig. 3 .Fig. 4 .

 34 Fig. 3. Top view of the experimental confi guration

Fig. 5 .Fig. 6 .Fig. 7 .

 567 Fig. 5. Example of time evolution of the rotation angle q

abFig. 9 .ZFig. 10 .

 910 Fig. 9. Surface grid of the blade. a M12e4 mesh, b M14e5 mesh

Fig. 11 .Fig. 12 .Fig. 13 .

 111213 Fig. 11. Horizontal slice through the M14e5 mesh. a Top view, b close up of the region above the blade

17

 17

Fig. 14 .Fig. 15 .Fig. 16 .

 141516 Fig. 14. Free surface at t = 2.7 s. a External view, b shaft view

Fig. 17 .Fig. 18 .

 1718 Fig. 17. Free surface at t = 3.9 s. a External view, b shaft view

Figure 21 Fig. 19 .Fig. 20 .Fig. 21 .

 21192021 Figure21plots the quantity d ip over half a blade length (denoted by L hb ). As point P is the localization where the resultant force acts, this fi gure shows that the migration of this point towards the shaft of the blade is accurately calculated. Finally, Figs. 22 and 23 concern the effi ciencies. At the beginning of the stroke, it can be seen that the main difference between h and h 0 is mainly

Fig. 22 . 5 qFig. 23 .

 22523 Fig. 22. Evolution of effi ciencies h f (a) and h p (b). The effi ciencies are defi ned in Eq. 5

Table 1 .

 1 Characteristics of the confi gurations studied

	Parameter	Sim1	Sim2	Sim3
	Mesh name	M12e4	M14e5	M14e5
	Cell number Time law	120 000 Adaptive	1 400 000 Uniform (dt = 0.005 s)	1 400 000 Adaptive
	Target Courant number	0.3	-	0.3

Acknowledgments. The authors gratefully acknowledge the scientifi c committee of CINES (project dmn2049) and IDRIS (project 0129) for the attribution of computer time.


Esubalew Alemneh
email: esubalew_jalew@etu.u-bourgogne.fr

Sidi-Mohammed Senouci
email: sidi-mohammed.senouci@u-bourgogne.fr

Philippe Brunet
email: philippe.brunet@u-bourgogne.fr

An Energy Efficient Smartphone Sensors' Data Fusion for High Rate Position Sampling Demands

Keywords: inertial navigation systems, GPS, traffic safety, pedestrian navigation, dead reckoning

Many smartphone-based traffic safety applications have been proposed in literatures. These applications demand very high position sampling to safeguard vulnerable road users. European Telecommunication Standards Institute (ETSI) has defined time interval between Cooperative Awareness Messages for collision risk warning to be between 1s and 0.1s. This implies that geographical awareness information has to be sampled between the frequencies 1Hz and 10Hz inclusive. However, an investigation we made depicts that current smartphones can't support such high rate location sampling. Even though they meet the aforementioned sampling requirements, high rate sampling of position data is an energy hungry process. In light of this, we have proposed an energy efficient position prediction method that fuses GPS and Inertial Navigation Systems (INS) sensors data to estimate pedestrians' positions at high rate. INS based dead reckoning is performed to extrapolate positions from last known location when GPS reading is unavailable and when GPS fix is realized the reading is used to correct dead reckoning parameters in addition to serving as a location fix. The proposed solution is compared to a position prediction method that relies solely on GPS data on two selected pedestrian trajectories. The result demonstrates that fusing GPS and INS position data has an average improvement of 30% and 61.4% in error in distance and direction respectively. The proposed position prediction algorithm is also applied to sensors data that are obtained by relaxing sampling rates with the objective of sparing smartphone's energy. In this regard, first energy efficiency of different position sampling rates of GPS and INS sensors are evaluated and then the algorithm is applied to the sampling frequencies that are proven to husband energy. The outcome of the evaluation is that the battery life of smartphones can be doubled by compromising accuracy of estimated distance and direction by only 11.5% on average.

I. INTRODUCTION

Smartphones have become part of everyday life of human beings. Gaming, navigation, monitoring environmental and traffic conditions, supervising health of possessors, intervention on our behalf, and social networking are some of their applications [START_REF] Xia | The power of smartphones[END_REF] [START_REF] Campbell | From Smart to Cognitive Phones[END_REF]. In recent days, due to their wide existence and increased storage and computation capability smartphones have become important ingredients of infrastructure-less traffic safety, infotainment and traffic efficiency applications. A staggering increase in their capability is further enhanced by their sensing potential which is another feature that makes the devices more intelligent and essential part of life.

Though the number of sensors may vary from one smartphone to the other, typical smartphone contains accelerometer, a digital compass, a gyroscope, a GPS, quad microphones, dual cameras, near-field communication, a barometer, light, proximity, and temperature sensors [START_REF] Campbell | From Smart to Cognitive Phones[END_REF]. However, smartphones are far from becoming apotheosis sensing devices. Producers of the devices purposefully make them compact and small so that they will be cheaper and fulfill user preferences. Size reduction and subjection of the devices to non-ideal environments compromise the sensor accuracy and reliability highly and this has become hindrance to fulfil application requirements. For instance, ETSI has enacted maximum and minimum time interval between beacon generations to be 1s and 0.1s for collision risk warning systems. Despite the fact that smartphone-based solutions are proposed for traffic safety by warning about-to-occur collisions based on positions of vehicles and vulnerable road users (VRU) [START_REF] Alemneh | PV-Alert: A fog-based architecture for safeguarding vulnerable road users[END_REF] [START_REF] Hernandez-Jayo | V-Alert: Description and Validation of a Vulnerable Road User Alert System in the Framework of a Smart City[END_REF], present smartphones' position sampling rate do not meet the minimum time interval positions have to be sampled.

A research conducted by Zhizhong et al. in 2013 [START_REF] Ma | Experimental evaluation of mobile phone sensors[END_REF] has indicated that maximum GPS sampling rates of smartphones is 1s [START_REF] Ma | Experimental evaluation of mobile phone sensors[END_REF]. Our investigation on maximum sampling rate of latest smartphone showed that things haven't changed till now. Losses in connections and imperfect behavior of the data collection application further increase the sparseness of the data [START_REF] Jaiswal | Location Prediction with Sparse GPS Data[END_REF]. From these facts we can conclude that location updates of smartphones do not meet applications requirements since pedestrians with maximum 1.83m/s walking and 4.2m/s running speeds [START_REF] Hamdane | Issues and challenges for pedestrian active safety systems based on real world accidents[END_REF] could face traffic accident in the time interval between two consecutive locations updates as vehicles impact velocity could reach up to 80m/s [START_REF] Han | Analysis of Vulnerable Road User Kinematics Before/During/After Vehicle Collisions Based on Video Records[END_REF]. Therefore, to celebrate smartphones' involvement for traffic safety we have proposed a solution that uses inertial sensors to predict positions in case of unavailability of GPS location fixes. The method integrates GPS and low-cost inertial sensors to satisfy high position sampling demand of traffic safety applications.

INS can be used to pinpoint position of smartphone holders and they are not dependent on any external measurements like radios [START_REF] Correa | A Review of Pedestrian Indoor Positioning Systems for Mass Market Applications[END_REF]. It is possible to extrapolate locations from reference points by calculating distance travelled and direction of movement from inertial sensors data. INS sensors data can be sampled at different rates. High sampling rates of GPS and INS sensors is known to drain smartphone batteries rapidly. Therefore, applications that demand fast position sampling must address the issue of energy consumption in order not exacerbate energy hunger of smartphones. This can be done by slacked sampling of GPS and INS sensors data instead of sampling them at their fastest possible time. However, we have to take a caution not to harm accuracy of position prediction by reducing the number of data points for estimation. Our work in this research includes evaluating energy efficiency of various GPS and INS sensors data sampling rates and application of the prediction algorithm proposed for sampling rates that are found to increase energy efficiency of smartphones. This paper reports an energy efficient dead reckoning method to estimate locations of VRUs particularly pedestrians. The locations are estimated at required rate from known location and kinematics of the pedestrian by integrating GPS and INS sensors data. The contributions of the research are;

x We have conducted an experiment to determine if different location methods endorsed by present smartphones fulfil high location sampling requirements of traffic safety applications. The experiment shows that maximum sampling rates of current smartphones' location methods do not meet application requirement of traffic safety applications.

x In order to cope with high location sampling rate demands of traffic safety applications, we have proposed a novel INS and GPS sensors data fusion method to track pedestrians using their smartphones. The method doesn't depend on map information or movement history of VRU. The evaluation of the solution using two selected trajectories shows that it is able to improve accuracy of position prediction that solely depends on GPS.

x Energy consumption of various GPS and INS sampling rates are experimentally investigated. The proposed position prediction method is applied to sampling rates which contributes to smartphones energy efficiency. The result shows that with small compromise of accuracy large amount of energy can be gained by relaxing sampling period of INS sensors.

The rest of the paper is organized as follows. Section 2 summarizes related works in INS based dead reckoning and energy efficiency. An investigation conducted to verify maximum sampling rate of location methods is presented in Section 3. Section 4 discourses detail of the proposed solution while Section 5 is dedicated to discussions of the evaluation results made on the proposed algorithm and energy efficiency of different GPS and INS data sampling rates. Finally, conclusion is presented in Section 6.

II. LITERATURE REVIEW

Traffic safety applications necessitate high rate position sampling to safeguard pedestrians and other VRUs from traffic accidents. Though smartphones have been recommended for traffic safety applications as sources of positions of VRUs [3-4] [10-11], their position sampling rate is limited. Therefore, at the time GPS fix is not possible position of VRUs must be dead reckoned using supplementary positioning systems. Dead reckoning (DR) can be applied in offline mode [START_REF] Kourogi | A method of pedestrian dead reckoning for smartphones using frequency domain analysis on patterns of acceleration and angular velocity[END_REF] or online mode [START_REF] Nawarathne | Dead reckoning on smartphones to reduce GPS usage[END_REF]. Offline DR is used for post-processing of locations after collecting known locations. Online DR is suitable for realtime systems where estimation has to be made immediately after a known location is captured. Most position estimation solutions proposed are either for vehicular navigation [START_REF] Nawarathne | Dead reckoning on smartphones to reduce GPS usage[END_REF][START_REF] Kao | Integration of GPS and dead-reckoning navigation systems[END_REF][START_REF] Liao | The Performance Analysis of Smartphone-Based Pedestrian Dead Reckoning and Wireless Locating Technology for Indoor Navigation Application[END_REF] or for indoor pedestrian navigation [12] [16-17]. In order to extrapolates next positions from known GPS locations for pedestrian traffic safety, an online DR that can be applied in indoor or outdoor environments is required.

Energy saving by lowering sensor data sampling, absence of GPS location signals, and a need to improve location accuracy are some of the reasons for applying position prediction methods in navigation. Outdoor navigation systems that are used in open sky scenarios can use GPS based locations.

However, to reduce smartphone energy consumption by lowering GPS sampling rate, GPS based localizations can be supplemented by dead reckoning that involves low energy cost sensors [12-13] [18]. INS sensors helps to improve positioning reliability in situations where GPS location is not available due to obstruction of GPS satellites by high buildings, trees etc. INS can be also integrated with GPS for best location accuracy [START_REF] Kao | Integration of GPS and dead-reckoning navigation systems[END_REF] as both methods suffer from their own drawbacks. According to Zengshan Tian et al. in [START_REF] Tian | Pedestrian dead reckoning for MARG navigation using a smartphone[END_REF] navigation information obtained from the INS sensors is boasted with continuity, high dataupdating rate, good short-term accuracy, and stability. In this research INS sensors are used to supplement GPS based positioning to fulfill higher location sampling rate requirements than the current smartphones can support. As our experiment has proved and other researchers have indicated [START_REF] Ma | Experimental evaluation of mobile phone sensors[END_REF] maximum GPS sampling rate is 1Hz which is below the minimum sampling frequency enacted by ETSI.

To make position predictions, distance travelled and direction of movement have to be estimated. Most dead reckoning solutions use accelerometer to predict displacement by first detecting step events and estimating stride length [16] [20-21]. Direction of movement is usually estimated from magnetometer and gyroscope data and position prediction is made after every stride or more. For problems that require very fast position predictions multiple estimations must be made in a period of one stride. Therefore, step detection and stride length estimation-based position prediction do not fit for such application. Hence, the proposed prediction method involves estimation of displacement by using conventional kinematic equations of distance and velocity. GPS readings are used to correct dead reckoning parameters periodically in addition to serving as correct positions of pedestrian trajectory. Accelerometer and magnetometer are used to estimate attitude of movement of pedestrian. Intermediate positions estimated from displacement and heading are set to pass through Kalman Filter (KF) for further improvement of the positions.

Gathering GPS and INS sensor data at high data-updating rates is an energy hungry process [START_REF] Ma | Experimental evaluation of mobile phone sensors[END_REF] [START_REF] Oshin | Improving the Energy-Efficiency of GPS Based Location Sensing Smartphone Applications[END_REF][START_REF] Huang | A power efficient pedestrian touring scheme based on sensor-assisted positioning and prioritized caching for smart mobile devices[END_REF][START_REF] Liu | Balancing Power Consumption and Data Analysis Accuracy Through Adjusting Sampling Rates: Seeking for the Optimal Configuration of Inertial Sensors for Power Wheelchair Users[END_REF][START_REF] Shih | Energy-Aware Pedestrian Trajectory System[END_REF]. [START_REF] Oshin | Improving the Energy-Efficiency of GPS Based Location Sensing Smartphone Applications[END_REF][START_REF] Huang | A power efficient pedestrian touring scheme based on sensor-assisted positioning and prioritized caching for smart mobile devices[END_REF] have shown sampling GPS at different sampling periods has different energy consumption. According to [START_REF] Ma | Experimental evaluation of mobile phone sensors[END_REF] [START_REF] Liu | Balancing Power Consumption and Data Analysis Accuracy Through Adjusting Sampling Rates: Seeking for the Optimal Configuration of Inertial Sensors for Power Wheelchair Users[END_REF] normal INS data sampling took 25% less battery power than fastest sampling. We have studied how much energy hungry different GPS and INS data samplings are at a very small sampling period. The application of position prediction method proposed on sensor data obtained by energy efficient sampling periods demonstrate that significant amount of energy can be saved with small effect on accuracy of prediction. The following section presents an experiment conducted to check maximum sampling rate of current smartphones.

III. MAXIMUM GPS SAMPLING RATE OF SMARTPHONES

Maximum GPS sampling rate of smartphones depend on factors like used location method, make of the phone, embedded GPS chipset, as well as version and type of operating system of the mobile device. For collision risk warning system ETSI enacted maximum and minimum time interval between beacon generations to 1s and 0.1s, respectively [START_REF]Intelligent Transport Systems (ITS); Vehicular Communications; Basic Set of Applications; Definitions[END_REF]. This implies that location sampling frequency should be between 1Hz and 10Hz since the beacons comprises of location information. To check if modern smartphones meet high GPS sampling requirements of collision risk warning systems we have conducted an experiment. The location methods supported by current smartphones are high accuracy, GPS only, and network (Wi-Fi a Mobile Internet). High accuracy location method, which is also termed Assisted-GPS (A-GPS), uses GPS chip and networks to get location fixes. GPS only location method uses solely GPS chips and network location methods may depend on Wi-Fi or mobile data internet connection.

All applications except our experimental application were denied location access and the mobile device was set to run only the default manufacturer processes. Samsung smartphone (Model-SM-A520F) was used exclusively for this controlled experiment. An Android application is developed which can read smartphone location at the highest possible rate. By enabling each location methods turn by turn the experiments were conducted in indoor and outdoor environments. Keeping the smartphone with the application running on front jacket pocket one of the researchers walked for about an hour for each experiment. The position data is saved in a log file. The average minimum possible sampling period is presented in TABLE 1. The highest average GPS fix frequency is achieved by GPS only location system. Using this location method, we can get location updates approximately each second. This value confirms with the result obtained in [START_REF] Ma | Experimental evaluation of mobile phone sensors[END_REF]. Location obtained from mobile internet and Wi-Fi based location methods' is too sparse for traffic safety applications. A-GPS sampling rate is lower than GPS only based location method as it uses network location methods when GPS fixes are not possible. The experiment reveals that GPS sampling rates of current smartphones do not satisfy application requirement of traffic collision risk warning systems. Thus, location methods that can provide positions at higher rates should work together with GPS based location to meet high rate location demand of traffic safety applications. The next section is dedicated to describe details of our solution that enables high frequency positioning.

IV. POSITION PREDICTION SYSTEM

Position prediction system proposed is aimed to extrapolate locations in a very short time interval (100ms). A GPS fix gained every second is used as correct location and to adjust dead reckoning parameters for forthcoming predictions. If no new location is received position estimation is made using INS sensors data. Components of the system are shown in Fig. 1. Each component is elaborated in subsequent sub-sections.

A. Sensors and Sensor Data Processing

Accelerometer and Magnetometer are two INS sensors used for the position predictions in addition to GPS sensor. Speed from GPS readings is used to estimate displacement by adjusting velocity of pedestrians obtained from accelerometer. Raw acceleration data is split into two by filtering out gravity. This process results in linear acceleration and gravity virtual sensors. The two software-based sensors are synthesized through a combination of accelerometer and gyroscope inputs and Android performs optimized filtering and sensor fusion processes by applying Extended KF [START_REF]Sensor Types[END_REF]. Magnetometer and gravity sensors are employed to predict movement direction of smartphone holder. As magnetometer is susceptible to working machines and metals in the vicinity to remove the highfrequency noise from its output, a low-pass filtering has been applied. After repeated experiments, we found that at cutoff frequency f = 0.01Hz the heading result is more durable and accurate in the environment the test is conducted.

B. Position Prediction

After pre-processing of sensor data, position is estimated by figuring out displacement and heading of the pedestrian [16] [START_REF] Tian | Pedestrian dead reckoning for MARG navigation using a smartphone[END_REF]. Linear acceleration is used together with speed obtained from GPS to estimate displacement. Heading is computed from accelerometer without motion component and corrected magnetometer values. Once the two values are available, intermediate position can be calculated using equation [START_REF] Xia | The power of smartphones[END_REF].

Where, (X 0, Y0) and (X, Y) are initial and intermediate positions, D is the displacement and β is the heading. The destination position is obtained by passing intermediate position throw KF to improve accuracy of prediction. KF equations presented in [START_REF] Alemneh | An Online Time Warping based Map Matching for Vulnerable Road Users' Safety[END_REF] are used since the method have improved GPS accuracy by 12.59%. The detail of calculation of displacement and heading are presented beneath.

C. Displacement Calculation

Calculating displacement using inertial based pedestrian dead reckoning usually involves step detection and step length estimation [17] [19]. This method is not applicable if position estimation has to be made in a very short time interval since many predictions may be needed in a period of one stride which is 0.4m to 1.0m in length [START_REF] Kang | SmartPDR: Smartphone-Based Pedestrian Dead Reckoning for Indoor Localization[END_REF]. This is true for the application in consideration as position prediction has to be made every 0.1s. Conventional kinematics equations shown on equation ( 2) are used to calculate displacement.

Where, Dt, Vt, and At are displacement, velocity and acceleration at time t, t is change in time, V0 is the initial velocity. V0 is set to 0 at the beginning. Velocity is calculated by integration of acceleration found from accelerometer without gravity. Since acceleration and velocity take vector values, integration is made for each axis. Due to accumulated error velocity extracted from linear acceleration tends to increase exponentially over a short period of time, see Fig. 2. Therefore, velocity from GPS fixes is fused and used as initial velocity at each GPS fix. Let V gps be speed obtained from GPS, and β be heading, then X and Y components of initial velocities at ground truth point are computed using equation (3) [START_REF] Kang | SmartPDR: Smartphone-Based Pedestrian Dead Reckoning for Indoor Localization[END_REF].

X = X0+Dsin(β) Y = Y0+Dcos(β) (1)
Vt = V0+At*t Dt = Vt*t + 0.5* At*t 2

(2)

We can't use only the velocity extracted from GPS fixes as it wouldn't capture frequent changes as GPS are sampled every one second and displacement is calculated every 0.1s. Fig. 2. Velocities from Accelerometer and GPS Readings and Fusion of the two. Integrated velocity is more stable in comparison with velocity from accelerometer. By taking this velocity the distance between two GPS fixes is found to be approximately the same as sum of distances between GPS points.

D. Heading Estimation

Heading is the deviation of the smartphones Y-axis from magnetic north, measured clockwise in the East-North plane. Corrected values from magnetometer and gravity portion of accelerometer are integrated to estimate attitude angle of the smartphone [START_REF] Zhou | Pedestrian dead reckoning on smartphones with varying walking speed[END_REF]. From the two set of values a rotation matrix is produced to map points from local coordinate system (LCS) to global coordinate system (GCS) to resolve random placement of smartphones. GCS is a real-world North-East-Gravity direction coordinate system. If Rt is the rotation matrix and γ is an azimuth which is degrees east of magnetic north in LCS, then the orientation is converted to GCS using equation ( 4).

The actual heading of a user is forecasted by accounting declination angle. Refer [START_REF] Kourogi | A method of pedestrian dead reckoning for smartphones using frequency domain analysis on patterns of acceleration and angular velocity[END_REF] [17] to get detail on heading estimation. Fig. 3 displays comparison of heading calculated using INS sensors and GPS bearing during a short trip in the test site. The test data which is collected from rectangular trajectory depicts that the two values are similar with only minor differences.

Comparison of Heading Estimation with Ground Truth GPS Bearing

V. EVALUATIONS AND RESULTS ANALYSIS

In this section evaluation result of the proposed method is presented, energy efficiency of various GPS and INS sampling rates are assessed, the proposed method is applied to energy efficient sampling methods and the results are discussed.

A. Evaluation of the Proposed Method

The proposed method is tested with a data collected using Android application that reads GPS points at the highest possible sampling period. The algorithm running on the smartphone calculates the displacement and heading angle and estimates locations at 0.1s interval. Accelerometer and magnetometer were set to sample at maximum rate and average values are taken for the estimation period. Ground truth data is collected using real-time kinematics (RTK) position device (AsteRx4 OEM). Two trajectories are conceived to test the proposed system. The first trajectory is a straight path with only one turn while the second is a rectangular path. The proposed GPS and INS Sensor based Prediction (GISP) is compared in terms of error in distance of predicted points and deviation of heading from ground truth points of walking path with Simple Linear Prediction (SLP) that uses only GPS data. The method utilizes speed and bearing of the most recent GPS fix to estimate the next position. Fig. 4 depicts average and maximum distance errors of predicted points from ground truth points for the two prediction methods on the two trajectories. Maximum prediction distance error is included to show the worst-case inaccuracies of prediction algorithms. For both trajectories the method that integrates GPS and INS sensors data has lesser errors than the one that uses only GPS data. In fact, integrating inertial sensors with GPS improved position prediction accuracy by about 30% percent since direction of movement is updated frequently from INS sensors. The prediction is more accurate for Trajectory 1 due to its simplicity. Fig. 5 shows average error in heading of the two methods. Error in the second trajectory is higher than the first due to more turns in the path. In GISP heading of the pedestrian at each prediction point is feed from INS sensors. Because of this, the method improves accuracy of heading estimation by about 61.4% from the one which relies only on GPS readings. Because of space limitation we can't include visualization of the prediction in this report.

B. Sampling Rates and Energy Consumption

Maximum sampling period of current smartphones are found to be one second. Higher sampling periods however could be used for energy saving as long as noticeable energy saving can be achieved and accurate prediction can be made. Energy consumption of GPS sampling at one, two, three and four seconds is evaluated. Relaxed sampling periods of INS sensor data can also help in saving energy. Two default Android and one user defined inertial sensor sampling rates are assessed. The details of inertial sensor sampling types considered are displayed in TABLE 2. Position prediction method explained in previous sections samples GPS and INS data at highest possible frequencies. Sampling frequency of GPS and INS affect energy consumption of smartphones and accuracy of position predictions [5] [22-24]. The higher sensor data sampling rate implies the higher energy consumption and the more accurate position estimation. In this section energy consumption of various GPS and INS sampling rates are assessed. Next, prediction exactitude of prediction method proposed using various sampling methods that are found to enable obtrusive energy saving are evaluated.

Three Android applications are developed to collect energy consumption information of the smartphone used for the experiment. Application1 is used to collect battery information of a phone that involve sampling neither GPS nor INS. Application2 collects battery detail of a smartphone that samples GPS at 1, 2, 3 and 4 seconds. Application3 is used to get battery information of a phone that samples GPS each second and INS at rates specified in TABLE 2. A smartphone with battery capacity of 3000mAh is used for the experiment. The following configurations are made before conducting the experiment. The smartphone is reset to factory setting. During testing no application runs except Application1, Application2 or Application3 and the smartphone is used strictly for this experiment. All types of connections are turned off and screen brightness is set to the lowest possible value. Root access of the android smartphone is obtained so as to read current discharge of the smartphone. The smartphone is charged to 100% before each experiment. For each experiment the smartphone holder walked for 50 minutes in outdoor and indoor walking paths. For additional 10 minutes the smartphone is placed on a flat desk in office.

The result of energy consumption of smartphones at different GPS Sampling period and INS Sampling types is indicated in Fig. 6. GPS reading at GPS sampling periods considered has no significant energy consumption difference. As stated in [START_REF] Oshin | Improving the Energy-Efficiency of GPS Based Location Sensing Smartphone Applications[END_REF][START_REF] Huang | A power efficient pedestrian touring scheme based on sensor-assisted positioning and prioritized caching for smart mobile devices[END_REF] location updates at lower intervals have very subtle differences in the battery decay over time. In onehour test, battery level of the smartphone has decreased from 100% to 97.17% in 1s GPS sampling period and from100% to 97.84% in 4s GPS sampling period. The drop of battery level due to GPS sampling in our experiment in general is modest since the smartphone is used solely for this purpose and it has high battery capacity. Fastest INS sampling consumes energy quicker than the other two sampling types. A battery that could have stayed for 23 hours with GAME and 26 hours with USER DEFINED INS sampling types will run out only in 11 hours.

C C. Sampling Rates and Prediction Accuracy

Small-scale energy overhead on different GPS data sampling period doesn't convince us to compromise accuracy by choosing higher GPS sampling periods. In case of INS data sampling on the other hand a noticeable energy saving can be achieved by selecting energy efficient INS data sampling type. Accuracy of position prediction depends on whether the prediction is made on straight or on curved trajectories. Hence, the evaluation is made for straight and curved sections of Trajectory 2. From the rectangular trajectory three curves and four straight lines are considered. In addition to GISP and SLP, predictions that involve GAME (viz GISP_GAME) and USER DEFINED (viz GISP_USER) inertial sensor sampling are compared.

Average error in distance and attitude angle of predicted positions for straight and curved trajectories is displayed in Fig. 7 and8, respectively. GISP is the most accurate prediction method on both straight and curved trajectories since it takes many INS data points to make the prediction. On straight line there is no abrupt change of direction and this is the reason why SLP performs well. In curved trajectories SLP performs worst because of lack of updated movement direction until the next GPS fix is attained. GISP_USER do not work well on straight and curved roads since the number of inertial sensors data taken on average are two and sometimes the prediction may be made using only one or no data points. Unless we take considerable number of sensor data, due to noisy nature of inertial sensors prediction is poor. Prediction accuracy of GISP_GAME is the second best to GISP. Sampling period of this prediction is lower than GISP but higher than GISP_USER.

VI. CONCLUSION

Even though smartphones have become integral part of traffic safety applications, they do not support high position sampling requirement of the applications. Therefore, we have proposed inertial navigation system-assisted positioning system that works together with GPS sensors to predict locations of pedestrians at a very high rate. The method uses GPS fixes to correct dead reckoning parameters and as ground truth points for the next position estimations. A method that uses only GPS readings to extrapolate positions is compared with the new method using data from two trajectories conceived. As the proposed system integrates inertial sensors that have very high short-term accuracy and GPS which has very good long-term accuracy it outperforms GPS only based prediction. So as to make the position prediction energy conserving, energy efficiency of various data sampling rates of GPS and INS data sampling are investigated. Results show that by sparsening GPS samples the energy that could be saved is minimal while accuracy could be highly affected. However, by reducing sampling frequency of INS data sampling obtrusive energy is able to be saved with minor accuracy compromise. Therefore, inertial sensors can be integrated with GPS sensor in an energy efficient way to satisfy high position sampling rate demand of traffic safety applications.

Fig. 1 .

 1 Fig. 1. Component of Proposed Position Prediction System

Fig. 3 .

 3 Fig.3.

Fig. 4 .

 4 Fig.4. Average and Maximum Error in Distance of Predicted Points for the two Trajectories.

Fig. 5 .

 5 Fig. 5. Error in Heading of GISP and SLP for the Two Trajectories

Fig. 6 .

 6 Fig.6. Energy Consumption of Different GPS and INS Sampling Rates

Fig. 7 .

 7 Fig.7. Average Error in Distance of Different Position Prediction Types

Fig. 8 .

 8 Fig.8. Average Error in Movement Direction of Pedestrian in Different Position Prediction Types GISP_GAME is recommended position prediction method when pedestrians are in non-traffic risky regions since it doubles the life of battery with 11.5% compromise of accuracy.

TABLE 1 .

 1 AVERAGE MINIMUM LOCATION SAMPLING PERIOD OF DIFFERNT LOCATION METHODS

TABLE 2 .

 2 ANDROID INS SAMPLING TYPES

		FASTEST	GAME	USER DEFIND
	Average Sampling	5	20	50
	Period (ms)
	Average Number	20	5	2
	of Samples