0 Introduction

With the rapid integration of unmanned aerial vehicles (UAVs) into modern military operations, autonomous aerial refueling (AAR) has attracted increasing interests from many research institutions and scientific laboratories. There are two typical AAR methods, the flying-boom refueling (FBR) and the probe-drogue refueling (PDR)^[1]. PDR is able to refuel multiple aircrafts simultaneously with low cost and simple way. In PDR, a tanker aircraft deploys a flexible hose behind it with a drogue on the end of the hose, and its dorking phase is affected by tanker wake atmospheric turbulence and receiver fore-body. PDR requires skillful pilots to maneuver a thin probe into a moving drogue with a small diameter at an acceptable closure rate. It has always been a challenging operation for both manned and unmanned missions because of mishaps and damages easily caused by any tiny deviation. It is essential to improve PDR process to enhance the safety of both tanker and receiver. Therefore, the research on dynamic modeling and active control plays an important part in the AAR technology^[2].

There are three approaches in PDR dynamic modeling: Partial differential equation method, finite element method based on mechanics of materials, and lumped parameter principle with rigid link kinematics. The former two are complex, with a large amount of calculation, which could hardly meet real-time requirement. Kamman et al^[3-5]. studied the modeling and simulation work of PDR for many years who have gained significant achievements. They presented a hose-drogue assembly model with a series of rigid links connected with frictionless ball-and-socket joints based on lumped mass, finite segment method. This model preferable reflected the dynamic behavior of hose-drogue system but the length of hose is constant. Hu et al^[6-7]. took the hose as a particle system, and calculated the hose balance position in terms of static equilibrium conditions, and studied the motion of hose-drogue influenced by atmospheric turbulence. Wanget al^[8]. proposed a variable length of multi-level ideal simple pendulum model based on Kamman′s work, and carried out comprehensive analysis of the features of hose-drogue system. Nevertheless, the methods above adopt newtonian mechanics method and need additional constraints to calculate internal force. Wei et al^[9-11]. studied the relationship between the drogue position and interference forces due to receiver fore-body based on a link-connected model and analyzed the drogue dynamic characteristics qualitatively based on system identification method, and proposed a lower order dynamic model to describe drogue dynamics under the bow wave effect.

Obviously, these studies, as well as the reports of flight tests, pointed out that the drogue position is susceptible to disturbances consisting of atmospheric turbulence, trailing vortex of the tanker and bow wave effect caused by the receiver. Among them, the receiver bow wave has a minimal effect on the drogue position compared to tanker wake and atmospheric turbulence^[12]. And the drogue motion is turned out to be instability especially in the atmospheric turbulence. The instability is able to cause an unpredictable drogue movement and increases the difficulty of a successful AAR operation. Hence, it is essential to stabilize the drogue in the presence of disturbance. The active control issues in the new context have been raised and initially investigated by some researchers. Kuk et al.^[13-14] developed an active drogue stabilization approach by installing a set of four foldable wings around the rear hub cap and performed comprehensive dynamic wind tunnel experiments with a scale model prototype manufacturing, also implemented a feedback control law to avionics hardware design. Williamson et al.^[15] presented a controllable drogue equipped with canopy materials attached to the end of the strut-arm basket to generate required side force. Moreover, a control law was executed to attenuate effects due to offsets from a desired position. Liu et al.^[16] proposed a boundary control scheme based on Lyapunov′s direct method to regulate the hose′s vibration and designed a disturbance observer to estimate the input disturbance. Wang et al^[17] designed an integral sliding mode backstepping controller for the hose whipping phenomenon relying on a permanent magnet synchronous motor and high-precision position sensors. Since all these concepts are feasible propositions lack of practical design details in the public domain, the control strategy for the hose-drogue system deserves our further exploration.

This paper is based on our previous work^[18] and is distinguished from it. We propose a hose-drogue dynamic model. The main differences between them are as follows: (1) The restoring moment is taken into account for more exact external force analysis; (2) a part of the formula of force expression is written in vector form for more explicit; (3) conceptual active control strategies using both PID and LQR methods for the hose-drogue system are investigated based on the hose-drogue dynamic model, and the control effects are well compared.

1 Method and System Description 1.1 Principle of Kane method

Kane equation method provides a unified formula to solve the problem of rigid multi-body dynamics without requiring solution of the constraints equations, which is better at formulating equations of motion for complex mechanical system than classical approaches such as Newton method. Kane equation method employs generalized velocities to depicted the system, which are a group of variables independent from each other and can be used to uniquely define the velocities and angular velocities of the system. The number of generalized velocities equals to that of independent dynamical equations, as well as the degrees of freedom (DoFs) of the system^[19-20].

For a system with n rigid bodies and m DoFs, Kane equation can simply be expressed as a set of scalar equations

$ {F_r} + F_r^ * = 0\;\;\;\;\;r = 1, \cdots ,m $

(1)

where, F_r and F_r^* represent the generalized active force and the generalized inertia force of the system with respect to generalized velocity u_r, respectively, which can be expressed as

$ \left\{ \begin{array}{l} {F_r} = \sum\limits_{k = 1}^n {\left( {{\mathit{\boldsymbol{F}}_k} \cdot {\mathit{\boldsymbol{V}}_{k,r}} + {\mathit{\boldsymbol{M}}_k} \cdot {\mathit{\boldsymbol{\omega }}_{k,r}}} \right)} \\ F_r^ * = \sum\limits_{k = 1}^n {\left( {\mathit{\boldsymbol{F}}_k^ * \cdot {\mathit{\boldsymbol{V}}_{k,r}} + \mathit{\boldsymbol{M}}_k^ * \cdot {\mathit{\boldsymbol{\omega }}_{k,r}}} \right)} \end{array} \right. $

(2)

where k stands for the index of rigid bodies. F_k(M_k) and F_k^* (M_k^*) represent the active force (moment) and the inertia force (moment) acting on the body k, respectively; V_{k, r}= $\partial $V_k/$\partial $u_r the partial velocity; and ω_{k, r}=$\partial $ω_k/$\partial $u_r, means the partial angular velocity of the body k with respect to the generalized velocity u_r; V_k the velocity of mass center of body k, and ω_k the angular velocity.

1.2 System description

In terms of lumped parameter method, the hose-drogue system is formulated under the following assumptions: (1) The hose is modeled as a series of constant length finite number of cylindrical shaped rigid bars connected with frictionless spherical joints; (2) the quality and external loads associated with each hose-bar are focused on the joints; (3) there is no rotation motion around central axis of hose due to the material characteristic of hose and exterior structure of bars; (4) the drogue is regarded as a mass point at the end of the hose. And the refueling pod moving with the tanker under the straight and level flight condition is connected with the other end of the hose.

The hose-drogue system dynamic model and coordinate systems are showed in Fig. 1, where S_g represents an inertial reference frame in which the system moves; S_p the trajectory frame of the tanker; S_k bar frame of each hose-bar. In Fig. 1, the origin O_k locates in the mass center of each hose-bar. X_k points to the tanker along the symmetry axis of bar; Y_k is orthogonal to symmetry plane and points right; Z_k points down in the symmetry plane. Based on the assumptions above and frames defined, there are two state arguments used to depict the relative position between each hose-bar and the tanker: Relative pitch angle θ_k1(the angle between the projection of O_kX_k in the plane of O_pX_pZ_p and O_pX_p in trajectory frame with downward as positive), and relative yaw angle θ_k2(the angle between O_kX_k and the plane of O_pX_pZ_p in trajectory frame with left deviation as position).

According to the angles defined and the coordinate systems, the coordinate transformation matrices among S_p and S_k can be obtained from

$ \mathit{\boldsymbol{S}}_p^k = \left[ {\begin{array}{*{20}{c}} {\cos {\theta _{k1}}\cos {\theta _{k2}}}&{\sin {\theta _{k2}}}&{ - \sin {\theta _{k1}}\cos {\theta _{k2}}}\\ { - \cos {\theta _{k1}}\sin {\theta _{k2}}}&{\cos {\theta _{k2}}}&{\sin {\theta _{k1}}\sin {\theta _{k2}}}\\ {\sin {\theta _{k1}}}&0&{\cos {\theta _{k1}}} \end{array}} \right] $

(3)

As suggested above, the hose-drogue system is assumed to be n level rigid bars connected to each other. Taking the typical hose-bar k for consideration, its spatial position can be depicted by θ_k1 and θ_k2. So the pose of the hose-drogue system can be expressed with number of 2n of independent generalized velocities expressed as $({{\dot \theta }_{11}}, {{\dot \theta }_{12}}, \ldots, {{\dot \theta }_{k1}}, {{\dot \theta }_{k2}}, \ldots, {{\dot \theta }_{n1}}, {{\dot \theta }_{n2}}) $.

2 Dynamic Model Formulation 2.1 Kinematics analysis

As shown in Fig. 1, the length of each hose-bar is assumed to be constant l_k, and p_k means the distance vector from hose-bar (k－1) pointing to hose-bar k. So p_k in S_p can be defined as

$ {\mathit{\boldsymbol{p}}_k} = {l_k}{\left[ { - {C_{k1}}{C_{k2}} - {S_{k2}}{S_{k1}}{C_{k2}}} \right]^{\rm{T}}} $

(4)

where C_ki=cosθ_ki, S_ki=sinθ_ki(i=1, 2). The position vector of hose-bar k in the trajectory frame can be expressed as

$ {\mathit{\boldsymbol{r}}_k} = {\mathit{\boldsymbol{r}}_{k - 1}} + {\mathit{\boldsymbol{p}}_k} $

(5)

As S_p and the refueling pod move, it is feasible to seek the position of all hose-bars by Eqs.(4), (5), respectively.

Then, the velocity and acceleration of hose-bar k in S_p can be found by the first and the second differentiating as

$ {\mathit{\boldsymbol{V}}_k} = {\mathit{\boldsymbol{V}}_{k - 1}} + {\mathit{\boldsymbol{\omega }}_k} \times {\mathit{\boldsymbol{p}}_k} $

(6)

$ {\mathit{\boldsymbol{a}}_k} = {\mathit{\boldsymbol{a}}_{k - 1}} + {\mathit{\boldsymbol{\alpha }}_k} \times {\mathit{\boldsymbol{p}}_k} + {\mathit{\boldsymbol{\omega }}_k} \times \left( {{\mathit{\boldsymbol{\omega }}_k} \times {\mathit{\boldsymbol{p}}_k}} \right) $

(7)

where ω_k and α_k are the angular velocity and angular acceleration of hose-bar k in S_p, defined as

$ {\mathit{\boldsymbol{\omega }}_k} = {\left[ {{S_{k1}}{{\dot \theta }_{k2}}{{\dot \theta }_{k1}}{C_{k1}}{{\dot \theta }_{k2}}} \right]^{\rm{T}}} $

(8)

$ {\mathit{\boldsymbol{\alpha }}_k} = {\left[ {{S_{k1}}{{\ddot \theta }_{k2}} + {C_{k1}}{{\dot \theta }_{k1}}{{\dot \theta }_{k2}}{{\ddot \theta }_{k1}}{C_{k1}}{{\dot \theta }_{k2}} - {S_{k1}}{{\dot \theta }_{k1}}{{\dot \theta }_{k2}}} \right]^{\rm{T}}} $

(9)

2.2 Parameters calculation for Kane equation 2.2.1 Partial velocity and partial angular velocity

Converting the velocity and angular velocity of hose-bar k into the form of generalized velocity in S_p as

$ \left\{ \begin{array}{l} {\mathit{\boldsymbol{V}}_k} = {\mathit{\boldsymbol{V}}_{k,j1}}{{\dot \theta }_{k1}} + {\mathit{\boldsymbol{V}}_{k,j2}}{{\dot \theta }_{k2}}\\ {\mathit{\boldsymbol{\omega }}_k} = {\mathit{\boldsymbol{\omega }}_{k,j1}}{{\dot \theta }_{k1}} + {\mathit{\boldsymbol{\omega }}_{k,j2}}{{\dot \theta }_{k2}} \end{array} \right. $

(10)

The partial velocity V_{k, j1}, V_{k, j2} and partial angular velocity ω_{k, j1,} ω_{k, j2} are calculated as follows

$ {\mathit{\boldsymbol{V}}_{k,j1}} = \frac{{\partial {\mathit{\boldsymbol{V}}_k}}}{{\partial {{\dot \theta }_{j1}}}} = \left\{ \begin{array}{l} {l_k}{\left[ {\begin{array}{*{20}{c}} {{S_{j1}}{C_{j2}}}&0&{{C_{j1}}{C_{j2}}} \end{array}} \right]^{\rm{T}}}\;\;\;\;\;\;\;j \le k\\ {\left[ {\begin{array}{*{20}{c}} 0&0&0 \end{array}} \right]^{\rm{T}}}\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;j > k \end{array} \right. $

(11)

$ \begin{array}{*{20}{c}} {{\mathit{\boldsymbol{V}}_{k,j2}} = \frac{{\partial {\mathit{\boldsymbol{V}}_k}}}{{\partial {{\dot \theta }_{j2}}}} = }\\ {\left\{ \begin{array}{l} {l_k}{\left[ {\begin{array}{*{20}{c}} {{C_{j1}}{S_{j2}}}&{ - {C_{j2}}}&{ - {S_{j1}}{S_{j2}}} \end{array}} \right]^{\rm{T}}}\;\;\;\;\;\;\;j \le k\\ {\left[ {\begin{array}{*{20}{c}} 0&0&0 \end{array}} \right]^{\rm{T}}}\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;j > k \end{array} \right.} \end{array} $

(12)

$ {\mathit{\boldsymbol{\omega }}_{k,j1}} = \frac{{\partial {\mathit{\boldsymbol{\omega }}_k}}}{{\partial {{\dot \theta }_{j1}}}} = \left\{ \begin{array}{l} {\left[ {\begin{array}{*{20}{c}} 0&1&0 \end{array}} \right]^{\rm{T}}}\;\;\;\;\;\;\;\;j = k\\ {\left[ {\begin{array}{*{20}{c}} 0&0&0 \end{array}} \right]^{\rm{T}}}\;\;\;\;\;\;\;\;j \ne k \end{array} \right. $

(13)

$ {\mathit{\boldsymbol{\omega }}_{k,j2}} = \frac{{\partial {\mathit{\boldsymbol{\omega }}_k}}}{{\partial {{\dot \theta }_{j2}}}} = \left\{ \begin{array}{l} {\left[ {\begin{array}{*{20}{c}} {{S_{j1}}}&0&{{C_{j1}}} \end{array}} \right]^{\rm{T}}}\;\;\;\;\;\;\;\;j = k\\ {\left[ {\begin{array}{*{20}{c}} 0&0&0 \end{array}} \right]^{\rm{T}}}\;\;\;\;\;\;\;\;\;\;\;\;\;\;j \ne k \end{array} \right. $

(14)

2.2.2 Inertia force and inertia moment

Based on the state vectors defined above, the inertia force and inertia moment can be calculated in frame S_p below.

$ \begin{array}{l} \mathit{\boldsymbol{F}}_k^ * = - {m_k}{\mathit{\boldsymbol{a}}_k} = \\ - \sum\limits_{j = 1}^k {{m_k}\left( {{\mathit{\boldsymbol{V}}_{k,j1}}{{\ddot \theta }_{j1}} + {\mathit{\boldsymbol{V}}_{k,j2}}{{\ddot \theta }_{j2}}} \right) - {m_k}{\mathit{\boldsymbol{\xi }}_k}} \end{array} $

(15)

where

$ {\mathit{\boldsymbol{\xi }}_k} = \sum\limits_{j = 1}^k {{l_k}\left[ {\begin{array}{*{20}{c}} { - 2{S_{j1}}{S_{j2}}{{\dot \theta }_{j1}}{{\dot \theta }_{j2}} + {C_{j1}}{C_{j2}}\left( {\dot \theta _{j1}^2 + \dot \theta _{j2}^2} \right)}\\ {{S_{j2}}\dot \theta _{j2}^2}\\ { - 2{C_{j1}}{S_{j2}}{{\dot \theta }_{j1}}{{\dot \theta }_{j2}} - {S_{j1}}{C_{j2}}\left( {\dot \theta _{j1}^2 + \dot \theta _{j2}^2} \right)} \end{array}} \right]} $

(16)

$ \begin{array}{*{20}{c}} {\mathit{\boldsymbol{M}}_k^ * = - {\mathit{\boldsymbol{J}}_k}{\mathit{\boldsymbol{\alpha }}_k} - {\mathit{\boldsymbol{\omega }}_k} \times \left( {{\mathit{\boldsymbol{J}}_k}{\mathit{\boldsymbol{\omega }}_k}} \right) = }\\ { - \left[ {{{\left( {\mathit{\boldsymbol{S}}_p^k} \right)}^{\rm{T}}}\mathit{\boldsymbol{J}}_k^k\mathit{\boldsymbol{S}}_p^k} \right]\left( {{\mathit{\boldsymbol{\omega }}_{k,k1}}{{\ddot \theta }_{k1}} + {\mathit{\boldsymbol{\omega }}_{k,k2}}{{\ddot \theta }_{k2}}} \right) - {\mathit{\boldsymbol{\eta }}_k}} \end{array} $

(17)

where

$ {\mathit{\boldsymbol{\eta }}_k} = {J_{ky}}\left[ {\begin{array}{*{20}{c}} {2{C_{k1}}S_{k2}^2{{\dot \theta }_{k1}}{{\dot \theta }_{k2}} + {S_{k1}}{S_{k2}}{C_{k2}}\dot \theta _{k1}^2}\\ { - 2{C_{k2}}{S_{k2}}{{\dot \theta }_{k1}}{{\dot \theta }_{k2}}}\\ { - 2{S_{k1}}S_{k2}^2{{\dot \theta }_{k1}}{{\dot \theta }_{k2}} + {C_{k1}}{C_{k2}}{S_{k2}}\dot \theta _{k1}^2} \end{array}} \right] $

(18)

In Eq.(17) J_k^k is the inertial matrix of hose-bar k in S_k. Assuming that the axial element J_kx of each hose-bar is zero due to the diameter d_k far less than the length l_k. Additionally, considering the axial symmetry of the configuration of the hose-bar, J_ky and J_kz are taken as the same values. The similar matrix form is used to denote the inertia matrix of the drogue since it is also axisymmetric. But its axial elements are finite value.

2.2.3 External force analysis

The resultant external force F_k acting on the center of gravity of hose-bar k, consisting of both gravitational m_kg, restoring moment due to hose bending R_k and aerodynamic force D_k, is

$ {\mathit{\boldsymbol{F}}_k} = {m_k}g + {\mathit{\boldsymbol{R}}_k} + \frac{{{\mathit{\boldsymbol{D}}_{k - 1}} + {\mathit{\boldsymbol{D}}_k}}}{2} $

(19)

A restoring moment within the hose that tends to straighten emerges when the hose is bent due to an applied forces, which can be modeled based on a simple beam theory. As Fig. 2 shows, a bent hose is represented by two adjacent hose-bars, and the hose restoring force can be treated as an equivalent external force R_k^[21].

The magnititude of the force R_k acting on hose-bar k can be written as

$ {\mathit{\boldsymbol{R}}_k} = \frac{{8EI\gamma }}{{{l^2}}} $

(20)

where I is the area moment of inertia of the hose′s cross section and E the modulus of elasticity; l the length of the hose-bars. They are calculated as

$ I = \frac{{{\rm{ \mathsf{ π} }}\left( {d_o^4 - d_i^4} \right)}}{{64}} $

(21)

$ \gamma = \arccos \left( {\frac{{{\mathit{\boldsymbol{p}}_k} \cdot {\mathit{\boldsymbol{p}}_{k + 1}}}}{{\left\| {{\mathit{\boldsymbol{p}}_k}} \right\|\left\| {{\mathit{\boldsymbol{p}}_{k + 1}}} \right\|}}} \right) $

(22)

where d_o(d_i) is the outside (inside) diameter of the hose.

The force of R_k is in the plane spanned by the hose-bar vectors p_k and p_k+1 with the direction of ${e_{{R_k}}} = \frac{{{\mathit{\boldsymbol{p}}_{k + 1}}}}{{{\mathit{\boldsymbol{p}}_{k + 1}}}} - \frac{{{\mathit{\boldsymbol{p}}_k}}}{{{\mathit{\boldsymbol{p}}_k}}} $.

D_k is composed of D_{t, k} and D_{n, k}that are tangential and normal aerodynamic drag separately, which are computed as^[4]

$ {\mathit{\boldsymbol{D}}_{t,k}} = \left\{ { - \frac{1}{2}\rho {{\left( {{\mathit{\boldsymbol{V}}_{t,k}} \cdot {\mathit{\boldsymbol{n}}_k}} \right)}^2}{\rm{ \mathsf{ π} }}{d_o}{l_k}{c_{t,k}}} \right\}{\mathit{\boldsymbol{n}}_k} $

(23)

$ \begin{array}{*{20}{c}} {{\mathit{\boldsymbol{D}}_{n,k}} = \left\{ { - \frac{1}{2}\rho {{\left\| {{\mathit{\boldsymbol{V}}_{t,k}} - \left( {{\mathit{\boldsymbol{V}}_{t,k}} \cdot {\mathit{\boldsymbol{n}}_k}} \right){\mathit{\boldsymbol{n}}_k}} \right\|}^2}{\rm{ \mathsf{ π} }}{d_o}{l_k}{c_{n,k}}} \right\}}\\ {\frac{{{\mathit{\boldsymbol{V}}_{t,k}} - \left( {{\mathit{\boldsymbol{V}}_{t,k}} \cdot {\mathit{\boldsymbol{n}}_k}} \right){\mathit{\boldsymbol{n}}_k}}}{{\left\| {{\mathit{\boldsymbol{V}}_{t,k}} - \left( {{\mathit{\boldsymbol{V}}_{t,k}} \cdot {\mathit{\boldsymbol{n}}_k}} \right){\mathit{\boldsymbol{n}}_k}} \right\|}}} \end{array} $

(24)

where n_k=[C_k1C_k2S_k2－S_k1C_k2]^T, ρ represents the air density; V_{t, k}=V_k－u_k, u_k the local air velocity which is dependent upon the tanker wake, steady flow as well as atmospheric turbulence. The air velocity relative to the tanker wake can be modeled using a Hallock-Burnham model^[22-23]; c_{t, k} and c_{n, k} represent the tangential and normal drag coefficients, which are involved with the Reynolds number of the local flow. Generally, they are obtained from experimental data.

The drogue is connected at the end of the hose, which is suffered the gravitational m_drog and aerodynamic force D_dro as well. So the resultant external force of the nth hose-bar is

$ {\mathit{\boldsymbol{F}}_n} = \left( {{m_n} + {m_{{\rm{dro}}}}} \right)\mathit{\boldsymbol{g}} + \frac{{{\mathit{\boldsymbol{D}}_n}}}{2} + {\mathit{\boldsymbol{D}}_{{\rm{dro}}}} $

(25)

The drag of the drogue can be computed as

$ {\mathit{\boldsymbol{D}}_{{\rm{dro}}}} = \left\{ { - \frac{1}{2}\rho {{\left( {{\mathit{\boldsymbol{V}}_n} - {\mathit{\boldsymbol{u}}_n}} \right)}^2}\left( {\frac{{{\rm{ \mathsf{ π} }}d_{{\rm{dro}}}^2}}{4}} \right){c_{{\rm{dro}}}}} \right\}\frac{{{\mathit{\boldsymbol{V}}_n} - {\mathit{\boldsymbol{u}}_n}}}{{\left\| {{\mathit{\boldsymbol{V}}_n} - {\mathit{\boldsymbol{u}}_n}} \right\|}} $

(26)

where d_dro is the diameter of drogue and c_dro the drag coefficient of the drogue which depends on the drogue′s geometric characteristics.

2.2.4 Dynamics equations of system

Substituting the partial velocity and partial angular velocity expressions along with the accelerations, angular velocities and angular accelerations obtained in Eqs.(11)—(26) into Eq.(2) gives a set of 2n number 2nd order nonlinear differential equations.

$ \sum\limits_{r = 1}^{2n} {{H_{pr}}{{\dot u}_r}} = {R_p}\;\;\;\;\;p = 1,2, \cdots ,2n $

(27)

where ${{\dot u}_r} $ are the derivative of generalized velocities, H_pr the functions of generalized coordinates, masses and inertia matrixes of each hose-bar and drogue, R_p the functions of generalized velocities and active forces, besides those of H_pr. Arranging these equations into a matrix form for computer simulation, the matrix [H_pr]_2n×2n is invertible, so that these equations are readily to be solved.

2.3 Validation of the model accuracy 2.3.1 Numerical solution procedure

Given that the physical and geometrical information of the hose-drogue system and motion parameters of the tanker as well as the initial angles of all hose-bars are known, the motion of the hose-drogue system may be tracked by the following steps.

(1) Define generalized coordinates and velocities, and set initial values;

(2) Compute the spatial position, velocity and angular velocity of each hose-bar;

(3) Compute the partial velocity and partial angular velocity with respect to each generalized velocity;

(4) Compute resultant extent force of each hose-bar;

(5) Establish Kane dynamic equations and solve with fourth-order fixed-step Runge-Kutta method;

(6) Compute the real-time position of drogue in inertial frame.

2.3.2 Steady-state characteristics

On account of the dynamic model put forward in the previous section, the simulation program is developed. The towing point of hose is assumed to be at the midspan of the tanker′s right wing, and the hose is apart into 20 segments. The values of simulation parameters needed are shown in Table 1.

In order to verify the correctness of the presented mathematical model, define drogue balance position as D=V_D/L_H, shown in Fig. 1. The steady-state drag and balance position at various flight conditions is calculated as shown in Fig. 3.

As Fig. 3 shows, the drogue drag and balance position are relative to the altitude and airspeed of the tanker. That is, the drogue distance (relative to the tanker) increases with the increasing tanker altitude, which is likely because the air density declines with increasing altitude and the drag of drogue decreases correspondingly. The steady-state position reduces while the tanker′s airspeed increases due to increased drag of drogue resulting from the speed up of tanker.

2.3.3 Dynamic characteristics

As reported in flight tests, atmospheric turbulence had great influence on the dynamic performance of the hose-drogue system. Even on a clear day, the drogue swings in the air with slight disturbance, hardly finishing refueling task. In this paper, the problem is studied by incorporating the Dryden wind-turbulence model^[24] into the hose-drogue model. Here, the turbulence generated is isotropic with the intensity of 1.5 m/s and a length scale of 530 m. The turbulence effects on the hose-drogue motion are shown in Fig. 4.

Fig. 4 shows the resulting responses, in which the system is suffered a 20 s duration of turbulence air velocity starting at t=10 s. In Fig. 4, as the tanker′s airspeed increases, the drogue′s displacement reduces. Since the speed of atmospheric turbulence remains still at fixed height, while the tanker′s airspeed is up, the relative rate of atmospheric turbulence′s speed with respect to tanker′s airspeed decreases. It appears that hose-drogue occurs irregular movement among the turbulence and it fails to settle to a steady-state position. It is consistently as described in the flight-test study.

3 Active Controller Design

Due to the influence by the atmospheric turbulence, the drogue moves in a continuous irregular pendulum way, which is difficult for receiver to dock accurately and safely. So it is essential to reduce the amplitude and frequency of drogue movement with active controller to improve the docking success rate. Here, a conceptual active control strategy by implementing aerodynamic control surfaces on the drogue for the hose-drogue system is investigated using the hose-drogue dynamic model.

3.1 Controller design model

The equations of motion for hose-drogue system derived in Section 2 are in the form of 2n nonlinear 2nd order differential equations. It is necessary to convert them into a set of 1st order differential equations so as to be solved numerically. And the Runge-Kutta 4th order method with a fixed-step size is used for this work. Thus, the 4n state variables of the set of angles and angular rates of the hose-bar orientation angles completely represent the current status. That is

$ \mathit{\boldsymbol{x}} = {\left[ {{\theta _{11}},{\theta _{12}}, \cdots ,{\theta _{n1}},{\theta _{n2}},{{\dot \theta }_{11}},{{\dot \theta }_{12}}, \cdots ,{{\dot \theta }_{n1}},{{\dot \theta }_{n2}}} \right]^{\rm{T}}} $

(28)

Additionally, the input vectors are external load terms including the horizontal force and the vertical force which are operated by the active control device. That is

$ \mathit{\boldsymbol{u}} = {\left[ {\begin{array}{*{20}{c}} {{u_{{\rm{hor}}}}}&{{u_{{\rm{ver}}}}} \end{array}} \right]^{\rm{T}}} $

(29)

In the matter of output variables, the horizontal and vertical positions of the drogue are of interest, which may be defined along with the state and input vector as

$ \mathit{\boldsymbol{y}} = {\left[ {\begin{array}{*{20}{c}} {{P_{ny}}}&{{P_{nz}}} \end{array}} \right]^{\rm{T}}} $

(30)

In order to carry out control law design for the system, the model above need to be linearized at an equilibrium point. Thus, the control design model in state-space format is obtained by numerically linearizing at a reference equilibrium condition as follows

$ \left\{ \begin{array}{l} {{\mathit{\boldsymbol{\dot x}}}_e} + \delta \mathit{\boldsymbol{\dot x}} = \\ \mathit{\boldsymbol{f}}\left( {{\mathit{\boldsymbol{x}}_e},{\mathit{\boldsymbol{u}}_e}} \right) + \left[ {\frac{{\partial \mathit{\boldsymbol{f}}}}{{\partial \mathit{\boldsymbol{x}}}}\left| {_{\mathit{\boldsymbol{x}} = {\mathit{\boldsymbol{x}}_e}\mathit{\boldsymbol{u}} = {\mathit{\boldsymbol{u}}_e}}} \right.} \right]\delta \mathit{\boldsymbol{x}} + \left[ {\frac{{\partial \mathit{\boldsymbol{f}}}}{{\partial \mathit{\boldsymbol{u}}}}\left| {_{\mathit{\boldsymbol{x}} = {\mathit{\boldsymbol{x}}_e}\mathit{\boldsymbol{u}} = {\mathit{\boldsymbol{u}}_e}}} \right.} \right]\delta \mathit{\boldsymbol{u}}\\ {\mathit{\boldsymbol{y}}_e} + \delta \mathit{\boldsymbol{y}} = \\ \mathit{\boldsymbol{g}}\left( {{\mathit{\boldsymbol{x}}_e},{\mathit{\boldsymbol{u}}_e}} \right) + \left[ {\frac{{\partial \mathit{\boldsymbol{g}}}}{{\partial \mathit{\boldsymbol{x}}}}\left| {_{\mathit{\boldsymbol{x}} = {\mathit{\boldsymbol{x}}_e}\mathit{\boldsymbol{u}} = {\mathit{\boldsymbol{u}}_e}}} \right.} \right]\delta \mathit{\boldsymbol{x}} + \left[ {\frac{{\partial \mathit{\boldsymbol{g}}}}{{\partial \mathit{\boldsymbol{u}}}}\left| {_{\mathit{\boldsymbol{x}} = {\mathit{\boldsymbol{x}}_e}\mathit{\boldsymbol{u}} = {\mathit{\boldsymbol{u}}_e}}} \right.} \right]\delta \mathit{\boldsymbol{u}} \end{array} \right. $

(31)

represented as

$ \left\{ \begin{array}{l} \mathit{\boldsymbol{\dot x}} = \mathit{\boldsymbol{Ax}} + \mathit{\boldsymbol{Bu}}\\ \mathit{\boldsymbol{y}} = \mathit{\boldsymbol{Cx}} + \mathit{\boldsymbol{Du}} \end{array} \right. $

(32)

Since the order of the state-space model is too large to compute and needs to be reduced. A reduction model with 6th or even 2nd order system is obtained by applying the balanced model reduction technique. The step response of linearized models are shown in Fig. 5. It is apparent that the lateral and the vertical motions of the drogue are decoupled from each other, thus a single-input- single output(SISO) design for each channel would be an option for active control design. The reduced order transfer functions of system are listed in Table 2. And the time domain behavior of full order and reduced order system are compared. It is without losing significant information in the system input-output behavior when the drogue moves less than 4.4 rad/s^[12]. The pole-zero maps of the linearized models are shown in Fig. 6. Fig. 7 shows the frequency response, and Fig. 8 the unit step response of the full order and reduced order models.

3.2 Regulator design

In order to devise a regulator that would maintain a relative stable position of drogue in the tanker wake and atmospheric turbulence, a control law is designed. A datum point is based on the natural steady-state position behind the tanker in the air, which is eithor set by the conclusion of steady-state characteristics analysis in Section 2.3 or estimated by mass of real flight tests and experiment datas. Assuming that the actual offsets of drogue position relative to the datum point can be measured. For instance, a differential global positioning system(DGPS) receiver or an electro-optical camera system on the drogue can be used. In order to keep the control design process simple, the 2nd order reduced model is chosen as the control design model. Furthermore, active controllers designed by proportion integral derivative(PID) and liner quadratic regulator(LQR) methods to compare the effectivities, respectively. The block diagram of the control simulation model is shown in Fig. 9.

3.2.1 PID regulator

The conventional SISO design techniques are using for lateral and vertical channels, respectively. The transfer functions of PID controller are explained as follows

$ {G_c}\left( s \right) = {K_{\rm{P}}} + \frac{{{K_{\rm{I}}}}}{s} + {K_{\rm{D}}} \cdot \frac{{Ns}}{{s + N}} $

(33)

where K_Py=10 787.14, K_Iy=10 476.72, K_Dy=1 854.79, N_y=105.76, and K_Pz=10 603.98, K_Iz=10 671.05, K_Dz=1 509.66, N_z=79.37.

3.2.2 LQR regulator

A LQR is also designed to track relative position using the state-space model, where integrated signal of outputs are added to augment the model so as to reduce the steady-state error, that is $\mathit{\boldsymbol{\dot z}} = y $. The augment model is as follows

$ \left\{ \begin{array}{l} \left[ {\begin{array}{*{20}{c}} {\mathit{\boldsymbol{\dot x}}}\\ {\mathit{\boldsymbol{\dot z}}} \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} {{\mathit{\boldsymbol{A}}_\mathit{\boldsymbol{r}}}}&0\\ \mathit{\boldsymbol{I}}&0 \end{array}} \right]\left[ {\begin{array}{*{20}{c}} \mathit{\boldsymbol{x}}\\ \mathit{\boldsymbol{z}} \end{array}} \right] + \left[ \begin{array}{l} {\mathit{\boldsymbol{B}}_\mathit{\boldsymbol{r}}}\\ 0 \end{array} \right]\mathit{\boldsymbol{u}}\\ \mathit{\boldsymbol{y}} = \left[ {\begin{array}{*{20}{c}} {{\mathit{\boldsymbol{C}}_\mathit{\boldsymbol{r}}}}&0\\ 0&\mathit{\boldsymbol{I}} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} \mathit{\boldsymbol{x}}\\ \mathit{\boldsymbol{z}} \end{array}} \right] \end{array} \right. $

(34)

where A_r, B_r, C_r are the 2nd order state-space model matrix. And the conventional cost function is used

$ J = \frac{1}{2}\int_0^\infty {\left[ {{\mathit{\boldsymbol{x}}^{\rm{T}}}\mathit{\boldsymbol{Q}}x + {\mathit{\boldsymbol{u}}^{\rm{T}}}\mathit{\boldsymbol{Ru}}} \right]{\rm{d}}t} $

(35)

where the matrix Q is a weight matrix relative to the error of state variables x, while R a weight matrix relative to the control effort.

According to the optimal control principle, the LQR controller is designed, and the gain matrix K is obtained by solving an algebraic Riccati equation as follows

$ {\mathit{\boldsymbol{u}}^ * } = - {\mathit{\boldsymbol{R}}^{ - 1}}{\mathit{\boldsymbol{B}}^{\rm{T}}}\mathit{\boldsymbol{P\hat x = }} - \mathit{\boldsymbol{K\hat x}} $

(36)

$ \mathit{\boldsymbol{PA}} + {\mathit{\boldsymbol{A}}^{\rm{T}}}\mathit{\boldsymbol{P}} - \mathit{\boldsymbol{PB}}{\mathit{\boldsymbol{R}}^{ - 1}}{\mathit{\boldsymbol{B}}^{\rm{T}}}\mathit{\boldsymbol{P}} + {\mathit{\boldsymbol{C}}^{\rm{T}}}\mathit{\boldsymbol{QC}} = 0 $

(37)

where Q=diag([12 000 8 000 8 000 10 000]),

R=diag0.40.5.

The closed-loop system of step input signal for PID and LQR controller are shown in Fig. 10.

The time response specifications for the hose-drogue system under consideration equipped with the proposed controllers are given in Table 3.

The simulation results, displayed in Fig. 10 and Table 3, clearly show that both controllers are successfully designed and PID controller exhibits better transient performance. Then the same PID controller law is also acted on the 6th order as well as the full order linear model, the similar simulation results are obtained as shown in Fig. 11.

3.3 Design evaluation

To further investigate the performance and validate the proposed controller, the control laws of PID and LQR are incorporated into the nonlinear model, respectively. The lateral and the vertical displacements of drogue are tested with and without the controller function, in which the system is suffered from the same turbulence air velocity as motioned in Section 2. As Fig. 12 shows, it is apparent that both PID controller and LQR controller greatly reduce the displacements. The maximum lateral and vertical displacements due to atmospheric turbulence reach 0.5—1 m with the controller powered off; while they are only 0.5—1 cm with controller enabled. These results well coincide with the current design requirements of keeping drogue motion within 10 cm. Besides, PID controller exhibits better control effect and higher control precision than LQR controller. From Figs. 12(e)—(f), it is clear that the effect of both control laws are positive. The drogue position trajectory is irregular in the atmospheric turbulence without controller; while a fixed point is obtained with the controller enabled, which greatly improves the safety and accuracy of docking. As Fig. 13 shows, the maxium actuation force of the lateral and the vertical are 70 N and 40 N for LQR controller, respectively; while they are 20 N and 30 N for PID controller, respectively. These should be practically realizable forces with equipments.

4 Conclusions

A dynamic model of a hose-drogue for aerial refueling system based on Kane equation is developed using the finite-segment and lumped-mass approach, which does not require solution of the additional constraints equations for internal force. Numerical simulation results verify the correctness of the model, and simulations of hose-drogue position relative to the tanker under different flight conditions are conducted to testify the steady-state characteristic. The dynamic characteristic suffered from moderate-level atmospheric turbulence is also studied. The simulation results show the hose-drogue is unstable in air in turbulence.

Active controllers of PID and LQR are separately designed to stabilize the position of drogue by obtaining a linear model for the hose-drogue system at a reference fight condition. Also, both control laws are implemented into the orignal nonlinear model for a closed-loop simulation test with atmospheric turbulence, which show both controllers are of positive effect to reduce the oscillation amplitude of lateral and vertical motions of drogue. And by comparison, the PID controller exhibits better response and higher control precision.

In this paper, the influence of receiver bow wave, the design of control wing structure and related measurement sensors are not taken into account. Some of these works are currently under consideration and others are need to be further investigated in the future.

Acknowledgements

This work was supported in part by the National Natural Science Foundation of China (No.61533008), the Fundamental Research Funds for the Central Universities (No.NZ2016104), and the Funding of Jiangsu Innovation Program for Graduate Education (No.KYLX15_0276).