Study of Underwater and Wave Gliders on the Basis of Simplified Mathematical Models

Both underwater and wave gliders are known as autonomous unmanned energy-saving vehicles which have recently found applications for monitoring the world ocean. The paper under consideration discusses simplified mathematical models of these platforms enabling the straightforward parametric investigation into relationships between their parameters and performance. In its first part the paper discusses equations describing the motion of an underwater glider (UG) in a vertical plane as a basis for derivations relating geometric, kinematic and hydrodynamic characteristics of UG and its lifting system with relative differential buoyancy and pitch angle. Obtained therewith are formulae for the estimation of the UG glide path speed, lift-to-drag ratio, range of navigation and endurance. The approach is exemplified for typical cases of the UG conceived as winged bodies of revolution and flying wings. The calculated results feature dependencies of the UG speed on its configuration and volume as well as on the angle of attack for different magnitudes of relative buoyancy. Also considered is an optimal mode of operation, based on the maximization of the lift-to-drag ratio. The second part of the paper is dedicated to the estimation of the thrust and speed of a wave glider (WG), comprising a surface module (float) and underwater module represented by a wing, with the use of a simplified mathematical modeling intended to clarify the influence of the parameters upon the performance of the WG. The derivations led to an equation of forced oscillations of the vehicle accounting for the interaction of the upper and lower modules, connected by a rigid umbilical. The exciting impact of progressive waves of a given length and amplitude is found through the calculation of the variation of a buoyancy force in accordance with the Froude–Krylov hypothesis. The derivatives of time-varying lift with respect to kinematic parameters, entering the equation of vertical motion of the WG, as well as coefficients of instantaneous and time-averaged thrust force, are found by resorting to the oscillating hydrofoil theory. The derivation of the available thrust and the approximate calculation of the drag of the vehicle with account of wave and viscous components enable the evaluation of the speed of the WG for the prescribed geometry of the craft and wave motion parameters.


Introduction
Both underwater and wave gliders are known as autonomous unmanned energysaving vehicles which have recently found applications for monitoring the world ocean. They can also be viewed correspondingly as underwater and interface agents of a global information robotic system.

Underwater Glider (UG)
The UG is an autonomous underwater vehicle propelled by gravitational forces produced through an excess of weight and/or volume, resulting in a positive or negative buoyancy. In terms of technology evolution, the UG is a further step from the so-called diving buoys, whereby, due to the installation of lifting elements, the vehicle acquires a Da Lyu et al. [15] explored the winglet effect on both hydrodynamics and the trajectory of a blended body UG. Xiacheng Wu et al. [16] investigated the influence of the wing position of the dynamics of UG.
Other endeavors regarding the UG configuration include the modeling and simulation of a disk underwater glider with issues covering its dynamic stability as well as a structural selection and numerical analysis of its pressure hull, Liu et al. [17], Zhao [18].
In Stryczniewicz et al. [19], the 3D CFD model was created to enable the simulation of a UG motion in six degrees of freedom for a small Reynolds number including the numerical integration of the equations of motion of an underwater glider shaped as a body of revolution with wings (BRW). In order to provide an adequate accuracy of the flow calculations, the authors independently meshed the domain surrounding lifting elements and the far-field around it, claiming that replacing the deforming meshes approach with the one using sliding meshes results in a high accuracy of force prediction.
Energy consumption estimation and minimization is one of the major issues discussed by the community of researchers implementing UG technology. Ming Yang et al. [20] applied the motion parameter optimization for different observation missions for Petrel-L UG in order to maximize the gliding range and minimize the energy consumption for one gliding cycle. They used the inner penalty function method and adopted a certain nonlinear variation rule of the buoyancy loss. The authors stated that the proposed gliding strategy could increase the gliding range by 24% for the UG Petrel-L.
Considered in Jiafeng Huang et al. [21] are a set of issues related to the design and motion simulation of a new UG for the depth of 400 m and a speed of 2 knots on the basis of a steady-state attitude formulation accounting for net buoyancy, and the displacement of a movable mass-block.
Shuxin Wang et al. [22] discusses the multi-disciplinary optimization of UG for improving endurance. The authors focus their attention on the design parameters which have significant influence on the gliding range for one cycle seen as the optimization target. These include the buoyancy factor, the compressibility of the pressure hull, hydrodynamic coefficients and motion parameters. The results, as claimed by the authors, show a significant increase of the gliding range. The analysis performed for the Petrel-L UG comprised matters of hydrodynamic configuration (body and wings geometry) and the pressure hull subsystem. The authors indicated that the Petrel-L UG designed with the proposed optimization approach was tested during two sea trials in the South China Sea.
It is worthwhile mentioning that the motion parameter optimization should account not only for the variation of net buoyancy but also for such effects as the deformation of the pressure hull and the variation of sea water density. Yang Y et al. [23] state that the buoyancy variation caused by these two factors is of the same order of magnitude as the nominal net buoyancy. Underlined in Yang Song et al. [24] is the fact that a full assessment of the energy loss due to the hull and water compressibility effects is quite important. Thus, the saved UG battery power can be used not only to feed on-board gauges, such as CTD, turbulence, acoustic and biochemical sensors (engaged in mesoscale eddy observation and turbulent microstructure and internal wave observation), but also some indispensable additional instrumentation, e.g., the Doppler velocity logger.
In Xue et al. [25], a formation of UGs has been considered which has advantages in sustained ocean observation with high resolution and the adaptation to complicated ocean missions. A multi-layer coordinate control strategy is developed for the fleet of hybrid underwater gliders to control the gliders' motion and formation geometry with optimized energy consumption. The feasibility of the coordinate control system and motion optimization method has been verified both by simulation and sea trials. The fleet of three Petrel-II gliders developed by Tianjin University has been deployed in the South China Sea, and the technology is available for larger groups of UG.

Wave Glider (WG)
The WG is an autonomous unmanned marine vehicle comprising two modules: a surface module (float), performing motions under the action of waves, and an underwater (winged) module connected with the latter through umbilical and comprising a wing system which generates thrust due to the motions of the upper module in the waves. The WG is a perfect interface platform designed to provide a gateway between the agents of the global robotic information system operating in the upper and the lower spaces. The platform is driven entirely by the renewable energy of the ocean, converting wave power into the translatory motion of the system and using solar panels fixed on the upper deck to feed instrumentation.
Wave-energy-harvesting principle has been implemented in marine robots exemplified by the wave glider (WG) of Liquid Robotics, introduced in 2005-2007, and an USV named AutoNaut, of the University of Southampton (2016). The Liquid Robotics-type WG has two modules. The upper module is subject to the action of waves and transmits its motions through a umbilical to the deeply submerged lower thrust-generating wing module which propels the vehicle, controlled by the rudder fitted on the surface platform. AutoNaut has a floating elongated displacement hull which is equipped with springed high-aspect ratio wings fitted at the extremities. Hereinafter we provide a short review of some relevant publications on the WGs.
Hine et al. [26] and Manley and Willcox [27] state that the key innovation of the Liquid Robotics Wave Glider (LRWG) is its ability to harvest the abundant energy in ocean waves to provide an essentially limitless propulsion. Ocean waves possess substantial power, and the WG harnesses this power to maintain an average reported forward speed of 0.78 m/s (1.5 kts) in seas with a 0.5-1 m wave height.
Note that the LRWG has been designed to withstand large open-ocean waves and strong winds with its low-profile surface float, high-strength tether and robust submerged glider. In 2007, during Hurricane Flossie, it demonstrated an ability to operate in seas with wave heights over 3 m and winds greater than 21 m/s (40 kts). While surveying along the Alaska coast, the LRWG successfully operated in waves over 6 m high and winds greater than 26 m/s (50 kts). The robustness of the LRWG design has led to many successful missions over long periods in strong and weak seas with many different deployment configurations.
Kraus and B. Bingham [28] defined their work as estimating the WG dynamics for precise positioning and developed a 2D model in longitudinal motion. They calculated the wave and driving force with the use of empirical data, and produced a simulation for some concrete values of wave height and period. Ryan N Smith et al.'s [29] study focuses on the Wave Glider Platform from Liquid Robotics and is targeted at determining a kinematic model for offline planning that provides an accurate estimation of the vehicle speed for a desired heading and set of environmental parameters. Given the significant wave height, ocean surface and subsurface currents, wind speed and direction, the authors present the formulation of a system identification to estimate the vehicle's speed over a range of possible directions. Based on the data collected over long-duration missions, it has been observed that the LRWG is able to maintain an average speed of approximately 0.8 m/s. Phillip Ngo et al. [30] applied Gaussian process models to existing wave data to predict the performance and carried out an effective method for forecasting WG velocity.
In Baoqiang Tian et al. [31], the nonlinear dynamic model of a wave-driven unmanned surface vehicle (WUSV) in two dimensions was established based on the analysis of its driving principle in the longitudinal direction. The authors calculated the wave and driving force, determined the hydrodynamic coefficients with the use of empirical data obtained with an experimental platform of WUSV, and, finally, presented the simulation results of the model.
To secure efficient wave-driven propulsion, one needs to use low-drag and stable upper displacing modules. Elhadad et al. [32] employed the Wigley model as the surface boat of the Wave Glider and calculated its resistance characteristics at a range of Froude numbers. Elhadad et al. [33] explored the overall stability performance of alternative hull forms of an Automated Oceanic Wave Surface Glider Robot using the Maxsurf software.
In Fuming Yang et al. [34] a numerical investigation of a Liquid Robotics type WG in head seas has been conducted with the use of a combination of commercial CFD packages Fine/Marine and STAR-CCM+. First, a simulation with the use of unsteady Reynolds Averaged Navier-Stokes (URANS built FINE/Marine with volume of fluid model-VOF) was executed regarding the flow around both the surface boat and the wing module. Secondly, using an overset mesh, a high-fidelity simulation of the passive eccentric rotation (PER) of the underwater hydrofoils was conducted. It was found that this rotation, coupled with a surge force acting upon the WG, is the main factor defining the propulsion efficiency of the vehicle. Not that the effect of the springs periodically returning wing elements to initial positions throughout oscillation cycles was modeled both in the experiment and in computations.
Feng et al. [35] established a cable model of a wave glider by considering connection characteristics such as a rigid rod, cable, multi-link chain and elastic rod. The results showed that the propulsion performance with different connection types of the wave glider was slightly different. Zhang and Xu [36] studied the motion relationship between the float and the sub glider by MATLAB.
Peng Wang et al. [37] studied a restricted circle-based position-keeping strategy for the wave glider. They found that under a given sea state, a smaller restricted circle ensures a better positioning accuracy, but requires more energy consumption. They also concluded that waves and currents both have significant impacts on the positioning accuracy and energy consumption for the wave glider position keeping.
Zhanfeng Qi et al. [38] numerically investigated the case of a semi-active flapping foil of the WG of NACA0012 cross-section with fully prescribed heaving motion and the pitching motion determined by the hydrodynamic force and torsion spring. Analyzed in this paper were the influences of reduced frequency, spring stiffness and critical pitching amplitude on the hydrodynamic characteristics. The authors found that propulsive performance of the foil depends on whether the reduced frequency is lower or higher than the so-called critical reduced frequency. In the former case, the propulsive performance of the flapping foil can be improved exponentially, whereas in the latter case the semi-active flapping foil cannot provide an effective thrust.
Andre Amador et al. [39] examined the use of four machine-learning frameworks to predict the speed response of the SV3 Wave Glider on timescales relevant to local motion planning. According to this paper, the models were trained using onboard measurements of incident waves and near-surface currents. The research draws on a 71-day dataset collected during two deployments conducted between October and December 2020 and April and June 2021 off the coast of Southern California. The observations span a wide range of environmental conditions, with significant wave heights and near-surface flow speeds ranging up to 6 m and 0.8 m/s respectively. The predictions generated by the bestperforming models showed good quantitative agreement when evaluated with multiple independent datasets, indicating that the selected features and models generalize well.
Rozhdestvensky and Zin Min Htet [40] developed a mathematical model of a ship with energy-saving wings, allowing to calculate the motions of the ship-plus-wings system induced by regular and irregular head waves. They also used the flapping wing theory to determine the wave-generated thrust of the wings, and, finally, to evaluate the energy efficiency design indices (EEDI) showing a reduction of carbon-oxide emissions due to the use of the renewable energy of the waves.
Kostas Belibassakis et al. [41] conducted a numerical and experimental investigation of an actively controlled dynamic bow wing for augmenting the ship propulsion in the head and quartering waves. Therewith, they produced a time-domain seakeeping analysis to estimate ship-foil responses.
The motivation for preparing this paper was to provide simple mathematical models for two most significant robotic agents of the global marine information system: UG and WG, which are propelled, correspondingly, by gravitational forces (net buoyancy) and the renewable energy of the ocean waves. Today there exists a lot of scientific literature on these vehicles, reporting on experimental, computational and full-size trial data. Recently, quite an effort has been applied to develop real-time simulation approaches which combine a full spectra of methods and environmental surveys. However, a simulation aiming to predict the functional behavior of the marine robots under discussion requires the use of appropriate transfer functions, i.e., the response of the system to all kinds of perturbations. Necessary for this purpose are analytical, asymptotic and computational methods of deterministic and probabilistic nature employing basic principles.
The paper includes: an introductory section containing a brief state-of-the-art review on the development of the UGs and WGs, a section on the UGs, including the problem formulation for a steady motion of an underwater glider with negligible variation of water density and no account for hull deformation under action of pressure, followed by an analysis of the influence of the design factors of the UG on its performance (speed on the trajectory and range). The WG section contains a formulation of the problem of motions of the two-module WG excited by regular waves, a solution of the problem and an estimation of wave generated thrust of the wing module. An approximate approach applied to the estimation of the speed of the WG is based on a hypothesis that uniform motion can be realized when the drag of the vehicle becomes equal to the available thrust for a given wave length.

Parametric Analysis of a Steady Motion of Underwater Glider in a Vertical Plane
For design and mission designation purposes, it is important to have a clear understanding of how geometrical, kinematical properties of the vehicle and its lifting system, as well as the relative excess of buoyancy, are related to speed on the trajectory and range for a given depth of operation.
Consider a steady mode of descent/ascent of an underwater glider (UG), equipped with a buoyancy engine (BE). Theoretically, such a mode can be realized when the projection of the excess of buoyancy force (positive or negative) onto the direction of the trajectory equals the drag force. For example, this is possible when the depth-wise gradient of density and the pressure hull compression are negligibly small, and for the case of compensational control of the buoyancy by means of the BE.
Consider a steady planing of a underwater glider in a vertical plane. Assuming the absence of acceleration, the equations of a dynamic equilibrium of the UG are written down, relating such parameters of ascent/descent as the vertical and horizontal components of the speed of the vehicle, angles of pitch and angle of trajectory with the configuration of the UG, position and characteristics of its lifting system, axial displacement of the battery and differential buoyancy volume.
The position of the UG during ascent/descent is shown in Figures 1 and 2. Therewith, the axes x and y of the earth's coordinate system are directed to the right and down, correspondingly. Additionally considered are: a body coordinate system x 1 0y 1 , whose axes are rotated with respect to x0y by a pitch angle β and a flow coordinate system x U 0y U . The following designations are introduced: α: angle of attack, γ = α + β: angle of trajectory, R W : weight force, R B : (Archimedes) buoyancy force, R x : drag force, R y : lift force, x p : abscissa of the point of application of the wing lift in the body coordinate system, (x 1W , y 1W ): coordinate of the center of gravity of the vehicle, x 1B : abscissa of the point of application of the buoyancy force, x 1Y : abscissa of the point of application of the lift force. Projecting the forces acting on the vehicle in steady motion mode onto the earth's system axes x and y, we obtain:  Horizontal axis projection Vertical axis projection Equilibrium equation for the moment with respect to the origin of the coordinate system.
It follows from Equation (1) that where ε(α) = R x /R y = 1/K: inverse hydrodynamic quality. It follows from Equation (5) that the angle of trajectory of the UG is inversely proportional to the arctangent of hydrodynamic quality of its lifting system. In other words, with an increase of the hydrodynamic quality of the UG's lifting system, its trajectory becomes flatter.
Accounting for the known relationships one can derive from (2) the following equation which shows that for the steady motion of the UG, its buoyancy (positive or negative) is countered by a resultant of the lift and drag force (see Figures 1 and 2). Write, alternatively, or On the other hand, it follows from the moment equation that Expressing R y in (11) with Equation (10), one can find where η represents the relative differential buoyancy, with η > 0 for the ascent and η < 0 for the descent. It is adopted above that the excess of buoyancy R B − R W = ∆R B = ρg∆V = ρgV 0 η is obtained by means of buoyancy engine (BE), which realizes a variation of buoyancy in the fore part of the UG for the same weight of the vehicle and that in the state of equilibrium, R W = ρgV 0 , where V 0 is the volumetric displacement of the UG in equilibrium mode on the surface.
With the Equation (10), one can derive an expression for the speed of the UG on the trajectory. Writing lift and drag through corresponding force coefficients where S wing indicates the reference area; taking into account that the inverse hydrodynamic quality can be written as ε(α) = R x /R y = C x /C y , we come to the following expression for the speed of UG: A ratio ηV 0 /S wing = ∆V 0 /S wing in (14) is associated with the volumetric loading upon the wing system. Its magnitude is significantly smaller for the UG of the flying wing (FW) type than for the UG configured as a winged body of revolution (BRW), which, similarly to the case of airplanes, results in a decrease of the speed on a trajectory when passing from fuselage configuration to flying wing configuration. Note that the drag coefficient entering (14) can be represented as a sum of the induced drag and viscous drag coefficients If the pitch angle equals β = 90 • (α = 0), the UG will execute a vertical ascent (descent). Therewith C y = 0 and C xi = 0, and the speed of vertical motion equals Compose ratio of the speed on trajectory and the speed of vertical displacement at or representing the induced drag coefficient in (15) in Prandtl format as where k i = C αα x /C α y 2 , and λ is the wing aspect ratio; rewrite (17) as The dependence of (19) as a function of lift coefficient C y for different magnitudes C x0 is shown in Figure 3, from which one can draw an important qualitative conclusion: for an UG with wings, its speed on the trajectory is less than the speed of vertical ascent/descent. Let us consider the relationship between the angle of trajectory and the angle of attack of the vehicle in more detail.

Rewriting the Equation (4) in the form
One can obtain the following quadratic equation with respect to the angle of attack: Figure 3 shows that the augmentation of the lifting capacity, characterized by the lift coefficient C y , entails a decrease of the UG speed on the trajectory.
The solution of the Equation (21) is straightforward: with the realized (real) magnitudes of the angle of attack α obtained for the variation of the angle of trajectory in the range in ascent mode, and in descent mode, where Note that Equation (12), derived above, interrelates the relative excess of buoyancy with the geometric parameters of the glider, as well as with its kinematic and hydrodynamic characteristics for realizable modes of steady gliding.
Evaluate the range of UG, assuming for simplification that the battery energy is mostly expended for the reversal of the propelling buoyancy force at a design depth H. It is easy to see that the energy necessary for a single-shot reversal of the thrust is equal to where, as earlier, V 0 is the volume of the pressure hull of the vehicle, ρ is the density of water, and η H is the efficiency of the pump. Here, the depth is introduced in m and the volume in m 3 . Therewith, the energy is obtained in joules. If the energy capacity of the UG batteries equals E b , then the density of energy (the energy contained in a volume unit of the battery) can be written as is the volume of the battery compartment of the vehicle. With a design depth H for the period of one cycle (single descent-ascent yo), the UG will clear in the horizontal direction a distance equal to 2Hctgγ = 2HK, where K = 1/ε represents the hydrodynamic quality of the UG. Then, with account of (1.26), one can approximately estimate the range of navigation of the UG using the formula where the nondimensional function f R = K/η characterizes the range of the UG. It follows from the expression of the range function f R that the range of the UG is directly proportional to its hydrodynamic quality K and inversely proportional to the relative buoyancy η. A joint consideration of the Equations (14) and (27) shows that the speed grows in proportion to the square root of the relative buoyancy η, whereas the range decreases in inverse proportion to η. If, for example, η is increased four times, then the speed on the trajectory would increase two-fold, whereas the range would decrease four times.
The time A to cover the distance R ug , associated with the endurance of the UG, can be found in the following way: where u 0 is the horizontal speed of displacement of the UG. Taking into account the Equations (4) and (14), as well as the expression for the range function f R , one can obtain after some simple derivations, where w υ = V 0 /S wing . Note that (29) gives the endurance A ug in seconds. It is of practical interest to consider a tuning of the vehicle to such a magnitude of the angle of attack which would result in the maximum hydrodynamic quality.
Following the classical theory of lifting surface, write the formula for hydrodynamic quality as Differentiating (30) with respect to the angle of attack and equating this derivative to zero, it is not difficult to find an angle of maximum hydrodynamic quality α opt , the corresponding (optimal) lift coefficient C yopt and the maximum hydrodynamic quality proper. We have If the optimal mode of gliding is realized, the corresponding speed on the trajectory, range and endurance can be found with help of the following formulae To exemplify the application of the approach introduced above, we consider configurations of the types of body of revolution with wings (BRW) and flying wings (FW) with the following basic parameters: BRW (ellipsoid of revolution of length l = 2.0 m and diameter d = 0.2 m with a wing of rectangular platform, span l w = 1 m, chord c w = 0.15 m and foil of 9% thickness, aspect ratio λ = 6.67), FW (triangular with sweep angle at the leading edge χ = 30 • , root chord c w = 0.705 m, span l w = 2.44 m and foil of 30% thickness, aspect ratio λ = 6.93). The initial volume of the pressure hull was assumed to be identical for both configurations and to equal V 0 = 0.042 m 3 . The derivatives of hydrodynamic coefficients associated with the lifting properties were taken from monograph [42], in particular: for BRW (C α y = 4.393, C αα x = 0.984) and for FW (C α y = 4.197, C αα x = 0.870). The coefficients C x0 of viscous drag at zero incidence were determined approximately: friction drag (by means of the method of equivalent flat plate and Scholz corrections for the thickness of the wing and the elongation of the ellipsoid of revolution) and pressure drag of the hull as a fraction of the drag of a disk of identical cross-section [43]. Reynolds numbers were calculated with the relevance of the corresponding characteristic length (length of the hull and the chord of the wing for BRW configuration, with the use of the chord of the wing for FW configuration). Friction drag coefficients were calculated depending on the regime of the flow (laminar or turbulent). Presented below are some results of the calculations. Shown in Figures 4 and 5 are graphs for speed on the trajectory versus the volume of the vehicle V 0 for different magnitudes of relative differential buoyancy η. Used therewith was a variant of the Equation (14) deployed in the following way:  of the UG for different values of relative net buoyancy for the same UG geometry. They show that for both the BRW and FW types of the UG even significantly enlarging the UG (self-similarly) does not result in an adequate increase of the speed. For the same volume of the pressure hull and relative net buoyancy, the UG of the FW type is slower on the trajectory. Figure 6 shows that both geometrical types of the UG augmentation of the angle of attack (i.e., the lift coefficient) lead to a slow-down of the vehicle on the trajectory for any given relative net buoyancy. This result is consistent with the general conclusion on the influence of lifting capacity on the UG speed, as illustrated in Figure 3.    Figure 8 is dependence of the function f R = K(η)/η, characterizing the range on the relative differential buoyancy. Therewith, for the UG of the FW type, we used (for a given η) a maximum value of hydrodynamic quality. The dashed line in Figure 8 corresponds to the maximum quality of the BRW type of the UG. However, as the angle of attack of the maximum hydrodynamic quality for the BRW type considerably exceeds the critical (separation) angle of attack, the range estimate shown, represented by a dashed line, is exaggerated. To approach a real situation in the case of BRW type, we show, with a solid line in the same graph, the result corresponding to the angle of attack close to its critical value (α = α sep = 11.4 • ). The general conclusion following from Figure 8: the ranges of navigation of the FW and BRW types relate to each other similarly to the ratio of their hydrodynamic qualities.  Figure 8 also illustrates the influence of the relative net buoyancy η (on the range function f R = K(η)/η). Based on the conclusions regarding the influence of the net buoyancy on the speed and range, one can state that for the UG the former changes in proportion to √ η and the latter-in proportion to 1/η. For example, increasing the net relative buoyancy four times would result in a two-fold increase of the speed, and, at the same time, would lead to a four-fold decrease of the range.

Shown in
To tune up a given UG for a given mode of functioning, one should know the relationship between the angle of attack α and the angle of trajectory γ = α + β. Such a function for configurations under discussion is plotted in Figure 9. One can easily see that FW-type gliders have a much flatter trajectory than the BRW-type gliders. At the same time, the practical angles of trajectory of the BRW type gliders fall into the range 20 • < γ < 40 • . The Figure 9 shows that the influence of relative net buoyancy η on the slope of the trajectory of both BRW and FW gliders is negligible. It also demonstrates that the smaller the angle of attack (lift coefficient), the flatter the trajectory. Another qualitative conclusion following from Figure 9 is that, other things being equal, the trajectory of the FW glider is much flatter than that of the BRW glider. In simple words, the BRW should be preferred for depth missions, and the FW should be preferred for range missions.
Thus, in Section 1, based on the main relationships of the steady gliding of UG in a vertical plane, the expressions were obtained to determine the speed of motion, the hydrodynamic quality, the range and the endurance, as well as the angle of trajectory for the UG of the BRW and FW types.

Estimate of Thrust and Speed of a Wave Glider on the Basis of a Simplified Mathematical Model
Proposed in the present paper is a method of estimation of thrust and speed of a WG based on a simplified mathematical model which employs the flapping wing theory [44][45][46]. Consider the vertical oscillations of a WG, with the upper module being represented with a float of a rectangular platform of length l and width b, both measured at a waterline plane. Let this float have a draft d and a waterplane area S 0 W L = lb with axes x and y of the right Cartesian coordinate system (x, y, z) lying in the plane of symmetry, whereas axes x and z lie in the plane of unperturbed water surface. Axis x is directed to the right, and axis y is directed upwards. In a general case of an inclined board of the float, the area of its cross-section by a plane parallel to x0z is described by a function S W L (y), therewith S W L (0) = S 0 W L . A sketch of two-module WG in waves in shown in Figure 10. Let progressive waves run to the right with a phase speed U phase = ω(k λ )/k λ , where in the general case of water of finite depth H, the time (circular) frequency ω(k λ ) is related to wave number k λ = 2π/λ (λ: wave length) by the following relationship, ensuing from a combined linear condition at the free surface of water If the WG is moving left with speed U 0 , the waves run relative to the float with a speed where Fr c = U 0 / √ gc is the Froude number based on the chord of the underwater module.
Therewith, the apparent frequency equals Assuming that the upper and lower models are rigidly connected, one can write down the equation of forced oscillations of the WG in the direction of axis y. Therewith, note that beside inertial forces, the WG is subject to a restoring time-dependent Archimedes force.
(ρ: water density, g: gravitational acceleration), unsteady lift force on the wing R ywing , due to its oscillations and unsteady exciting, Archimedes force R ywave , due to the periodic variation of the wetted volume of the float. The final equation has been obtained in the form: (M + m) ..
where M is the mass of the WG, m is the added mass of the non-lifting part of the WG, S W L (0) = dS W L /dy for y = 0.
Representing the lift force acting on the oscillating wing as an expansion with respect to kinematic parameters, and calculating the exciting force due to waves with the help of the Krylov-Froude hypothesis as a corresponding variation of the Archimedes force, where a w is the wave amplitude, write the WG forced oscillations equation as follows: To reduce Equation (45) to a non-dimensional form, introduce the half-chord of the wing c/2 as a typical length and speed U 0 of the WG as a typical speed. Then, the nondimensional time τ and the nondimensional vertical displacement of the float η can be written as follows: Introducing the instantaneous magnitude C ywing of the lift coefficient of oscillating underwater foil and its derivatives with respect to nondimensional kinematic parameters . η and ..

η,
(where S wing is the reference area of the wing), and also noting that on calm water the mass of the WG is counterbalanced by buoyancy force, that is M = ρS 0 W L d, where d is the draft of WG, one can rewrite the equation of the oscillations of WG in waves as follows: ..
k a can be considered as a Strouhal number based on the apparent circular frequency ω a .
To simplify, consider first the case of a straight-walled board float, for which It is easy to see that in this case the Equation (49) A particular (nonhomogeneous) solution η nh (τ) of Equation (49), satisfying the righthand side (δ w = 0), has the form As the homogeneous solution of (49), corresponding to free oscillations with damping decays exponentially, to calculate the thrust due to force oscillations of WG it is sufficient to consider nonhomogeneous solution (49). Therewith, the amplitude A w of forced oscillations has been found in the form The consideration of (57) reveals the possibility of resonant oscillations of WG, which are accompanied by a considerable increase of thrust.
With the aim of a maximum simplification of the mathematical model and the illustration of the proposed approach, we restrict ourselves to considering the case of purely vertical oscillations of the wing, whereby the thrust on the wing is completely defined by a so-called suction force. In this case, according to [31], the derivatives of the unsteady lift coefficient with respect to kinematic parameters . η and .. η can be calculated with the help of the formulae C . η y wing = −2πF(k a ), Functions F(k a ) and G(k a ), entering the Equations (58) and (59), represent, correspondingly, the real and imaginary parts of the Theodorsen function (see [44]), which can be expressed through Hankel functions of the second kind of the first order H 0 (k a ), in the following way: The coefficient c T of the thrust averaged across the period of oscillations and ideal efficiency of the wing-propulsor can be determined with the help of the formulae derived in [45,46]: Note that in (61) the amplitude of oscillations is related to the full chord of the wing, so that A * w = A w /2. The application of the derived simplified mathematical model enables one to interrelate the characteristics of waves and the WG proper with the thrust, generated on its underwater wing module, and, accordingly, with the speed of its motion.
To find the dimensional thrust averaged over the period of oscillations, one should multiply the thrust coefficient by the dynamic head of the flow and wing reference area, that is, Of interest from a design viewpoint is the estimation of speed U 0 of WG for given wave properties. To find U 0 , the thrust generated by the lower (wing) module should be equated to the drag of the WG. Noting that for steady motion this drag is a sum of the wave resistance of the float and viscous drag of the whole vehicle, one can write In this paper, for the purpose of an approximate evaluation, the wave resistance of the float is determined with the help of the formulae given in [47] for a shallow hull with wedge-like extremities, and viscous drag is found as a sum of the friction drag of the modules, calculated by the method of equivalent plate with account of the regime of the flow with addition of the drag of the umbilical cord.
For the calculation of the wave drag, take the following formula, relating a complex function of wave amplitudes A(θ) with wave resistance [48]: For the case of a straight-board hull of overall length l with wedge-like extremities and insert of constant width b and length l p at constant draft d, [47] gives the following expression for the calculation of the wave amplitude complex function A(θ): where i = √ −1, ν w = g/U 2 0 , C p = b/(l − l p ). Viscous drag has been estimated using the formula where Ω ra f t и Ω wing are wetted surfaces of the float and wing, correspondingly, C f ra f t and C f wing are friction coefficients, determined through formulae of Prandtl-Schlichting for the turbulent flow and Blasius formulae for the laminar flow In the Equation (68), the Reynolds numbers are composed using the length of the float and the chord of the wing, correspondingly. R x umb represents the drag of the umbilical cord connecting the upper and lower modules.
To validate the approach, a number of calculations have been carried out for a concrete variant of WG with the following parameters: l = 2 m, b = 0.6 m, d = 0.2 m, l w = 1.5 m, c = 0.3 m, m = 0.3 M, h = 5 m (h: submergence of the wing), ρ = 1000 kg/m 3 , ν = 1.3 × 10 −6 m 2 /s, d umb = h = 5 m, C x umb = 1 (assuming a laminar separation), and the wave length varied in the limits 3 m ≤ λ ≤ 50 m. Based on statistics of sea waves, the wave amplitude was assumed as a w = λ/40. In the calculation of the wave resistance, the length of the insert of constant width was adopted as l p = 0.97l.
Typical results of the calculation of curves of the drag and available thrust are presented in Figures 11 and 12, separately, for the wavelengths in the range λ = 10 m-50 m and in the range λ = 3 m-4 m, in order to demonstrate possible thrust resonant modes which for the input vehicle (float) length of 2 m are found to occur at wavelength λ = 3-see the thrust maximum around the UG's speed U 0 ≈ 0.34 m/s-and at λ ≈ 4 m-see the thrust maximum around the UG's speed U 0 ≈ 0.72 m/s. Naturally, for other lengths of the float, the resonant wavelength would be different.  The calculated data on thrust, speed of motion and ideal efficiency of the WG under study are summed up in Table 1. Both Figure 12 and Table 1 confirm that for short waves, there may occur a resonant growth of the thrust of WG.

Conclusions
In this paper the simplified mathematical models are discussed for the two important agents of the marine robotic information system, namely, for underwater gliders (UG) and wave gliders (WG). The former model is based on the equations of forces and moments equilibrium of a UG, moving (descending or ascending) along its trajectory in a vertical plane on the assumptions of depth-wise constant water density and zero compression of the pressure hull. The latter model reduces the problem of a WG as a two-module marine robot, performing motions under the action of regular waves, to finding oscillations of mass-spring type simple mechanical system excited by a harmonic force, and modeling the conversion of these wave-driven motion into the thrust of the lower (wing) module, propelling the vehicle forward.
Both models provide straightforward insights into the influence of the geometry and kinematics of the corresponding vehicle upon its performance and suggest guidelines for efficient design and more elaborate investigations based on computational methods and experimentation. For example, the UG model readily yields dependencies of speed on the trajectory and range upon the pitch angle, the geometry of the pressure hull and lifting system, the volume of the buoyancy engine, the relative volume of the battery compartment and the specific power of the batteries. The WG model enables one to evaluate the influence of wave parameters, dimensions and shape of the float and wing module on the efficiency of the wave-induced thrust and speed of the vehicle. It also allows for the investigation of options for tuning up the WG-plus-waves system for a chosen ocean basin.
The insufficiencies of the proposed approaches are due to excessive simplifications aimed at obtaining a better understanding of the functioning of the technology under study.
Further research envisages: • extending the proposed UG model to the case of depth-wise density variation and pressure hull compression with a more detailed description of its geometry and the use of computational mechanics; • extending the proposed WG model to the simulation of a multi-wing lower module connected to the float by means of a flexible and elastic umbilical.

Conflicts of Interest:
The authors declare no conflict of interest.