Precise Underwater Gliders Pitch Control with the Presence of the Pycnocline

Glider-based mobile currents observations are gaining increasing research attention. However, the quality of such observations is directly related to the pitch accuracy of the glider. As a buoyancy-driven robot, the glider will be strongly disturbed during the passage through the pycnocline. The pycnocline refers to the oceanic phenomenon where the density of the seawater changes abruptly with respect to depth. The presence of the pycnocline influences the pitch of the glider and consequently affects the quality of the observed currents data. In this work, we propose an actuator constrained active disturbance rejection controller (ACADRC) to improve the accuracy of the pitch angle control when gliders move across the pycnocline. For this purpose, the dynamical model of the glider is first derived. Then, the longitudinal plane motion model of the glider considering the density variation is analyzed. Based on that, we discuss three typical types of pycnocline encountered during glider profiling, which are the pycnocline, the inverted pycnocline and the multiple pycnocline. To alleviate the low accuracy of bang-bang control and proportion integration differentiation control, and furthermore, to mitigate the disturbance of pitch by sudden density changes, we propose the actuator constrained active disturbance rejection controller in conjunction with specific glider pitch actuator hardware constraints. Simulation results show that the proposed method has significant improvement in pitch control accuracy over the comparison methods.


Introduction
The underwater glider (UG) is a kind of buoyancy-driven oceanic observation robot [1][2][3]. They are not equipped with propellers, but rely on manipulating buoyancy and attitude to navigate [4,5]. Combined with the necessary GPS positioning and satellite communication, gliders are able to achieve autonomous observations for months. With low vibration, long range and high autonomy, UGs are widely used in oceanic sampling missions. Among them, glider-based currents observation has been brought into focus [6][7][8]. This is because the currents sampling via gliders greatly reduces the involvement of measurement vessels and collects data with high autonomy and acceptable quality. However, such missions are extremely demanding for pitch angle controls [9]. Prevailing UG current measurements utilize the acoustic Doppler to current profiler (AD2CP) sensor developed from Nortek [10]. The AD2CP has strict requirements on the beam angle, and the nonconforming attitude of the glider will bring obstacles to the subsequent data processing or even seriously affect the quality of the data. Therefore, to ensure the high precision control of pitch angle becomes the key issue for glider currents measurement. However, the following encountered challenges should be considered.
The sudden change of density in the pycnocline will bring significant disturbance to the pitch control of the glider [11]. Pycnocline is a common oceanic phenomenon that describes a stratified seawater structure [12]. Following the classical three-layer simplified model, the partial derivative of seawater density with respect to depth in pycnocline is significantly larger than that in the upper and lower layers. This results in a steep increase in seawater density over a limited depth range. The effect of pycnocline disturbance is particularly dramatic for the buoyancy-driven underwater gliders [13,14]. The deviation of the pitch angle from the expected value fluctuates significantly during the passage through the layer. Furthermore, the dynamics of gliders in the pycnocline needs further discussions.
The actuator constraints of the glider imposes challenges for the precise control of the pitch angle [15,16]. The adjustment of the pitch angle depends on the value of the buoyancy engine, the displacement of the movable mass and the hydrodynamic damping [17]. The buoyancy-driven system of the glider consists of the inner and the outer bladders. The adjustment of the buoyancy is achieved by the transfer of the oil between the two structures [18]. Usually, the buoyancy engine is set to a constant action strategy in order to dive to a well-defined depth. This means that the buoyancy engine is not controlled in real time during the pitch control procedure, but the density disturbance from the environment will be introduced in the dynamics modeling via the buoyancy engine [19][20][21]. In addition, due to the influence of hydrodynamics, the pitch of the glider will cause oscillation during the adjustment process, resulting in the decrease of control accuracy. So the displacement of the movable mass plays a key role in controlling the pitch of the glider. The position of the movable mass is driven by a drive system consisting of a DC motor and a gear train [22]. Therefore, the actuator constraints must be considered when designing the control inputs, because the actuator constraints limit the feasibility of the control inputs.
Another challenge for the precise pitch control of gliders is to achieve the suppression of disturbances by the control algorithm. Many outstanding researches have been carried out in the field of pitch control of underwater gliders [23][24][25]. Among them, the bangbang control has the concise structure. It does not depend on the control model [26]. For these reason, bang-bang control acts the basic method for glider pitch control. However, the controlled quantity of bang-bang is difficult to coordinate between overshoot and high accuracy. Proportion integration differentiation (PID) is widely used in glider pitch control for its engineering feasibility [27]. However, PID does not perform convincingly in the face of pycnocline control. In our previous work, applying the ADRC to control the glider pitch angle was introduced [22]. However, the previous work does not sufficiently consider the actuator constraints. Based on the above discussion, we believe the anti-disturbance algorithm for glider traversing the pycnocline considering the specific hardware demanding investigation.
In this work, we propose the actuator constrained active disturbance rejection method for controlling the pitch of the glider in the presence of the pycnocline. Main contributions are concluded as follows. Firstly, we derive the dynamics model of the glider. On this basis, we analyze the longitudinal plane nonlinear motion model containing the variable buoyancy term. Secondly, we propose the ADRC control algorithm considering actuator constraints to alleviate the interference of the pycnocline and meet the practical situation of glider hardware. Finally, we propose the metrics for controlling pitch angle in the presence of the pycnocline. The simulation results show that the proposed method has significant improvement in control accuracy of pitch angle compared with the comparison methods.
The paper is organized as follows: Section 2 derives the glider dynamic model. Section 3 demonstrates the ACADRC. Section 4 introduces three typical structures of the pycnocline and the control metrics. Section 5 presents the numerical simulation results and discussions. Finally, conclusions are summarized.

Dynamic Model Derivation
The derivation of the glider dynamics model has been presented in many outstanding works, such as [17,28]. This section introduces our up to date work compared to [3]. We simplified the form of the glider mass matrix, thus facilitating the decoupling process of the matrix model. Moreover, the model is modified in terms of the correction of the buoyancy term and the accurate calculation of the added mass to better describe the effect of the buoyancy change caused by the pycnocline.

Model Derivations
The main parts of our developed underwater glider OUC-III are shown in Figure 1. To facilitate illustration, the inertial coordinate system, the body coordinate system, and the velocity coordinate system are introduced. Applying the Newton's law in the inertial system, the general form of the forces exerted on the underwater glider is obtained, which is shown in Equation (1). The hydrodynamic forces and moments of the glider are described in the current coordinate system, while the hydrodynamic forces and moments can be represented in the body system by the transformation matrix between the current coordinate system and the machine system. The gravitational forces and moments represented in the inertial system can also be converted to the body system to represent the external forces exerted on the glider in the same coordinate system. Then, through transforming the expressions into the inertial system, the force analysis of the glider is completed, as demonstrated in Equation (2). Among them, the r v can be calculated through The calculation of components of the forces and moments in the body frame are depicted in Equation (3), where the transformation matrix R T BI and R BC are described by Equations (4) and (5).
According to the Lagrangian modeling theory, we discuss the kinetic energy of the glider. The toral kinetic energy consists of the portion of the glider and the portion of the added mass (Equation (6)). The derivation of the total kinetic energy with respect to the generalized quantities q g yields the linear momentum and the angular momentum expressed in the body system (Equation (7)). Denote the generalized mass matrix as I M , then the derivatives of P and Π with respect to time can be obtained as shown in Equation (8). Equation (11) shows the six degrees of freedom of accelerations of the glider obtained by transforming Equation (8). The derivatives of linear momentum and angular momentum with respect to time is expressed with the help of the transformation relationship of the linear momentum and the angular momentum between the body coordinate system and the inertial system, see Equations (10) and (11). Moreover, the specific expression of theİ M is shown in Equation (12). Figure 1. The OUC-III underwater glider systems. As indicated by the serial numbers, they illustrate 1 the mechanical shell, 2 the sampling loads, 3 horizontal wings, 4 the CTD profiler, 5 the vertical wing, 6 the antenna, 7 the electronic control unit, 8 the actuating system of the movable mass and 9 the buoyancy engine, respectively.
Q = dP/dt dΠ/dt =İ M q g +q g I M V Up to this point, all parameters in Equation (9) can be calculated. Combined with the specific kinematic equations (Equation (13)) and the expressions of hydrodynamic angles (Equation (14)), the simulation of the glider motion can be realized using the Runge-kutta method.

Longitudinal Model Considering Varying Seawater Density
This section introduces the longitudinal plane model of the glider. Since the movable mass of the glider keeps neutral and non-rotate, we denote the r v has the form as r v = [r v1 , 0, r v3 ] T . Moreover, for the vertical motion, V 2 , p and r maintains zero. Thus, the cross term of the added mass and the inertial hydrodynamic term driven by these motions is zero.
The adjusted added mass cross term C f , leaves only the (2nd, 3rd) element c f 23 that is not zero. On this basis, we denote the M, J, C has the form as Equations (15) and (16) shows, M can be further expressed in Equation (16), where the elements are demonstrated in Equation (17). Taking the above results into Equation (9) and considering the characteristics of the glider's vertical motion, the 3D of glider dynamic Equation is obtained.
Then, the model for the vertical motion of the glider is described in Equation (20). This work considers the effect of pycnocline on glider pitch control. Thus, the density distribution of the seawater cannot be considered as uniform. We corrected m b to reflect the abrupt change of density, refer to Equation (21). In Equation (21), the ρ u represents the uniformly distributed density, while ρ v depicts the varying density. Furthermore, based on the partial derivative of density with respect to depth, we modifiedṁ b as Equation (22) shows, where dep denotes the depth.ṁ b and m b represent the buoyancy term in the uniformly distributed density field.

Parameters Determination
In our previous work [3], we presented the determination of hydrodynamic parameters for the OUC-III glider, and the results are listed in Table 1. According to the research in [17], parameters in Table 1 is determined by the mechanical shell. Once these parameters are calculated via the curve fitting techniques, the hydrodynamic forces and momentum are able to determine combined with the specific hydrodynamic angles and the velocities. However, the added mass appeared in Equation (20) changes with the density. Other physical parameters of the OUC-III glider are listed in Table 2. Table 1. Hydrodynamic parameters of OUC-III.

Parameters
Values In this work, we calculated the added mass of OUC-III using the hess-smith panel method [29]. Firstly, the glider surface is divided in a continuous triangular panels using finite element software. Then, the node information is resolved from the result of the division, and the area, center and normal of the surface elements are calculated based on the obtained nodes, as shown in Figure 2. Then, construct the linear system of equations about the intensity of the Rankine source σ, as Equation (23). Where V f is the value of the velocity at infinity, n represent the direction of calculation, the number of the system of equations is the number of elements. The resulting σ is brought into Equation (24) to calculate the velocity potential Φ. The additional mass can be calculated from Equation (25). Note that the added mass is related to the density of the fluid. So the motion of the glider in the pycnocline needs to consider the change of the added mass.
Equation (26) shows the added mass terms for ρ = 1000 kg/m 3 , and these matrices need to be corrected for fluid density during the glider profile.

The Actuator Constrained Active Disturbance Rejection Control
This section introduces the concept of the actuator constraints and the basic structure of ACADRC. This part also introduces the basic settings of the control methods involved in the comparison in Section 5.

The Actuator Constrains
As discussed before, the buoyancy engine is prone to interfered by the environment. Furthermore, the OUC-III glider is not equipped with the elevator, which makes the effective maneuverability of hydrodynamics not accessible. So, the key parameter for pitch control is the longitudinal displacement of the movable mass. The physical system of the glider's movable mass is shown in Figure 3. The planetary gear is driven by a DC motor fixed to the movable mass block. The axial displacement of the movable mass is adjusted by meshing with the screw rod. Assuming that the movable mass of the glider has a displacement r m1 along the axial direction, resulting in the overall center of gravity of the glider moving the distance of l = m m r m1 /m v . Considering that the metacentric height of the glider is h, then the steady value of the glider pitch angle manipulated by the movable mass can be expressed by Equation (27). In order to ensure the driving effect of the movable mass at the ultimate angle, the axial tension of the screw corresponding to the extreme pitch angle should be calculated. The calculation is shown in Equation (28), where the θ sm represent the extreme pitch angle. Accordingly, the desired torch of the screw under such case can be calculated referring to the lead angle of screw and the equivalent coefficient of friction as Equation (29) illustrated.
F r = m m g cos θ sm (28) Denote the reduction gear ratio of planetary gear as i, the maximum continuous output torque as T c , and the transmission efficiency as η. Then the continuous output torque of DC motor can be calculated as Equation (30). The type of DC motor is selected according to the T m , and the actuator constraint for glider pitch control can be calculated referring to the rated speed of the DC motor. The calculated v L is 1.3 × 10 −3 m/s, which means that the control algorithm of the pitch must follow the constraints of the calculated v L . The control algorithm cannot exert influence on the pitch angle according to the ideal control input, which affects the control effect, but is close to the actual physical system. Moreover, the movement of the movable mass relies on a linear potentiometer to decode the control signal. This also contributes to the actuator constraints.

The ACADRC
The active disturbance rejection control mainly consists of the tracking differentiator (TD), the nonlinear state error feedback (NLSEF) and the extended state observer (ESO) as shown in Figure 4. The TD section arranges the transition process for the reference input signal, and this design avoids the drastic tracking of the abruptly changing reference signal. At the same time TD generates the differential signal of the reference output. The formulas of the TD are show in Equation (32), where r 0 denotes the speed factor, h is the filter factor, and fhan represents the optimal control function. The ESO estimates the system output and the total disturbance in real time, and generates the system error by making a difference with the differential signal generated by the TD. By using the system output signal θ and the input signal u to track and estimate the system state and disturbances, it is the key to achieve the suppression of environmental disturbances. Its principle is listed in Equation (33). β 01 , β 01 and β 01 are the adjustable parameters and z 1 , z 2 , z 3 are the estimations of the pitch angle and its derivatives of each order, respectively. fal is a continuous power series function with linear segments near the origin [22].

ACADRC Controller
Actuator Constraints NLSEF can be expressed as Equation (34). The control output of the classical ADRC can be derived via the output u 0 combined with the estimation of the disturbance z 3 . However, this control output does not fully consider the actuator constraints. By mapping the output u to the action space of the movable mass of the glider and combining it with the properties of the potentiometer-based control of the actuator, the output of ACADRC can be obtained.

The Comparison Methods
In this work, to illustrate the control effect of the proposed ACADRC, we chose bangbang control and PID control as the comparison methods. They are both commonly used methods in underwater glider engineering applications.The idea of bang-bang control is to output the actuator adjustment signal by comparing the the current controlled state with the target. In OUC-III, the movable mass has the coordinates [−0.65, 0, 0] T meters in the body-frame. Furthermore, the adjustable range is −0.65 ± 0.07. In order to prevent the actuator from excessive oscillation under bang-bang control, we selected the range of ±1 degree from the target pitch as the termination band of control. The PID constructs the actuator target signal based on the error between the current pitch angle and the target pitch angle. The incremental PID control is selected in this work, and the principle is shown in Equation (35), where ∆θ(k) = θ(k) − θ(k − 1). The selected parameters are k p = 0.0005, k i = 0.0000015, k d = 0.0000008.

Discussion on Pycnocline and Control Metrics
This section introduces the structure of the pycnocline used in Section 5. We also introduce the metrics that depict the effect of pitch control in the presence of the pycnocline.

Three Typical Types of Pycnocline
Three typical structures of pycnocline are adopted in this work, which are the typicallydiscussed pycnocline with three stratified layer [30], the inverted pycnocline [31] and the dual-pycnocline (double pycnocline) [32]. To facilitate referring, denote them as pycnocline 1, pycnocline 2, pycnocline 3. Their density-depth curves are shown in Figure 5. The red dot in the figure depicts the boundaries of pycnocline. In the typicallydiscussed pycnocline, the partial derivatives of the density with respect to depth are smaller in both the upper and lower layers of the jump layer, which causes the density of seawater in these two layers to increase slowly with depth. However, in the jump layer, the density increases steeply within a limited depth, which leads to undesired changes in pitch angle control. Generally, the density of seawater increases with depth, but in the inverted pycnocline, the density of seawater underwater is less than that near the surface. This may lead to a steep drop in depth for the underwater vehicle and also cause interference to the pitch angle. The dual-pycnocline can be originally considered as a continuous typically-discussed pycnocline, whose density exerts disturbances on the glider in two depth intervals. It causes multiple oscillations of the pitch angle control. The specific parameters of the three pycnoclines used in this section are listed in Table 3, where the ub represents the upper boundaries, and lb illustrates the lower boundaries. The intensity demonstrates the partial of density with respect to the depth among the jump layer, while the ui and li depicts the intensity of the upper and lower layer of the pycnocline, respectively. Taking the pycnocline 1 as an example, we demonstrate the interference of the pycnoncline exerted on to the glider pitch. As Figure 6 shows, the green band illustrates the pycnocline. It can be observed that the pycnocline influences the pitch of the UGs significantly. Table 3. Physical parameters of OUC-III.

The Control Metrics of Pitch Control
As shown in Figure 7, the profile of the glider can be divided into four processes. They are the accelerating diving phase, the stable diving phase, the dive-surfacing switching phase and the stable surfacing phase. Among them, the accelerating diving stage is more aggressive in pitch angle control in order to reach the target pitch angle as soon as possible. After reaching flare point (1), the glider is controlled according to the diving target angle. At (3), the glider dives to the set depth and starts to adjust buoyancy and attitude to prepare for surfacing. After the adjustment of the switching process, the glider reaches the flare point (2), which means that the glider starts to control the floating pitch angle. The determination of the flare time is defined in Equation (36). The δ depicts the error to index. We also define the tm to illustrate the time when the glider reach the set depth, which can be found in Equation (37). To determine flare point (2) (denotes as t f 1 ), the −19.1 and the +δ are chosen. Together with the t m , the interval for continuous controlling the pitch angle during the dive phase can be selected. As shown in Equation (38). The mean and variance of the glider pitch angle in this selected interval can reflect the control effect of the specific control algorithm.  [1] (36)

Numerical Simulation Results and Discussion
To examine the ACADRC control effect, we conduct the motion simulation of the glider in the presence of three types of pycnocline. The bang-bang control and the PID control are introduced to compare the control effect. The diving depth of the glider is set to 200 m. Because the interference of the pycnocline on the pitch control is mainly focused, we consider the movable mass of the glider as non-rotatable. Therefore, the roll angle and the heading angle of the glider are all zero degree during the simulation. According to the simulation results. Although bang-bang control is a common control techniques for engineering applications. However, it leads to drastic oscillations in pitch angle and insufficient suppression in the face of disturbances caused by density jumps. However, the pitch angle under bang-bang control did not show significant overshoot. PID control has a significant improvement in accuracy compared to bang-bang control. However, the control effect performs unstably between the diving and floating. Furthermore, the control results show that the pitch angle has interval steady-state errors and is seriously affected by density disturbances. In addition, the pitch angle of PID control has significant overshoot. The pitch angle controlled by ACADRC is more stable than the comparison method. The pitch angle can be controlled around the target pitch angle. Furthermore, it has a better suppression effect on density interference than the comparison method. The disadvantage lies in that the ACADRC still has overshoot. However, the overshoot of both ACADRC and PID occurs in the surfacing phase and the glider does not reach the flare point. As defined in Section 4.2, the control metrics for the pitch control during the selected interval are listed in Table 4. The results indicate that the ACADRC has the comprehensive performance over the comparison methods in the aspect of lower variance and closer mean value to the target pitch. To better visualize the error of the controlled pitch with respect to the target pitch, we provide the Figures 11-13. For each control algorithm, we show error plots for the diving control phase and the floating control phase. According to the figures, the pitch angle under bang-bang control fluctuates up and down along the target pitch, resulting in poorer accuracy. In contrast, PID has high control accuracy in some stages, and even performs best in all three algorithms. However, the PID cannot achieve the globally accurate pitch control. The performance of ACADRC is weaker than PID in some areas, but the control effect is smoother in the whole control cycle and the suppression of density variation is better than the comparison method. It is worth noting that the glider cannot suppress the deviation of the pitch angle even if the movable mass is adjusted to the limit position due to the abrupt change of density in the pycnocline. At this point, all three discussed methods are unable to control the pitch angle to maintain the set angle. However, ACADRC can make a continuous estimation of the disturbance to compensate the control signal and make the pitch angle return to the set state quickly. In contrast, the PID control still suffers from oscillations and delays. Therefore, we believe that ACADRC is capable of accurately controlling the pitch in the presence of the pycnocline.

Conclusions
In this work, we have discussed the problem of accurate control of the pitch angle of the glider in the presence of the pycnocline. We have established the longitudinal model of the glider considering density variation based on the derived six-degree-of-freedom glider dynamics equations. On this basis, the influence of three typical pycnocline structures on the pitch of the glider has been considered. The simulation results have shown that the abrupt density change caused by pycnocline will bring significant disturbance to the pitch angle. The actuator constraints of the glider in conjunction with the hardware characteristics of the glider have been discussed. Under this constraint, the input of the control algorithm cannot be fully mapped to the actuator action, which has deteriorated the effectiveness of the control algorithm. The ACADRC has been proposed for this purpose, which, on the one hand, has allowed estimating disturbances in real time during the control process. At the same time, the output has been improved to better adapt to the actuator constraints for the purpose of achieving better pitch control results. The simulations have shown the advantages of the proposed method in terms of pitch control accuracy compared to bang-bang control and PID control.
Future work will focus on further suppression of pitch angle errors caused by pycnocline, as well as field tests.

Data Availability Statement:
The experimental data supporting this study are provided within this paper.