• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    Modeling and simulation of a time-varying inertia aircraft in aerial refueling

    2016-11-23 08:05:11WangHaitaoDongXinminXueJianpingLiuJiaolongWangJian
    CHINESE JOURNAL OF AERONAUTICS 2016年2期

    Wang Haitao,Dong Xinmin,Xue Jianping,Liu Jiaolong,Wang Jian

    Aeronautics and Astronautics Engineering College,Air Force Engineering University,Xi’an 710038,China

    Modeling and simulation of a time-varying inertia aircraft in aerial refueling

    Wang Haitao,Dong Xinmin*,Xue Jianping,Liu Jiaolong,Wang Jian

    Aeronautics and Astronautics Engineering College,Air Force Engineering University,Xi’an 710038,China

    Aerial refueling;Equations of motion;Hose–drogue assembly;Time-varying inertia;Whipping phenomenon

    Studied in this paper is dynamic modeling and simulation application of the receiver aircraft with the time-varying mass and inertia property in an integrated simulation environment which includes two other significant factors,i.e.,a hose–drogue assembly dynamic model with the variable-length property and the wind effect due to the tanker’s trailing vortices.By extending equations of motion of a fixed weight aircraft derived by Lewis et al.,a new set of equations of motion for a receiver in aerial refueling is derived.The equations include the time-varying mass and inertia property due to fuel transfer and the fuel consumption by engines,and the fuel tanks have a rectangle shape rather than a mass point.They are derived in terms of the translational and rotational position and velocity of the receiver with respect to an inertial reference frame.A linear quadratic regulator(LQR)controller is designed based on a group of linearized equations under the initial receiver mass condition.The equations of motion of the receiver with a LQR controller are implemented in the integrated simulation environment for autonomous approaching and station-keeping of the receiver in simulations.

    1.Introduction

    Air to air refueling is an effective method of increasing the endurance and range of aircraft.Most recently,there has been increasing interest in autonomous aerial refueling(AAR)for the continuing research into unmanned aerial systems.Over the last decade,there have been a wealth of research and academic publications on the theoretical and practical aspects of automating the refueling process covering aircraft control,sensor systems,and their integration.1–4In order to further develop and evaluate these researches,it is critical to model the whole aerial refueling system with sufficient accuracy,taking into account all major factors including the aircraft,aerodynamic and atmospheric disturbances,and refueling apparatus.For the probe-drogue refueling,three of the most significant factors are the dynamic model of the receiver aircraft with the time-varying mass and inertia property,a hose–drogueassembly(HDA)dynamicmodelwiththe variable-length property and the wind effect due to the tanker’s trailing vortices.

    It is generally accepted that it is sufficient to consider the dynamics of an aircraft at a number of fixed,specified,overallweights.5If a mass change is within about 5%of the beginning mass after a period of 60 s,the constant mass assumption is considered acceptable.6Therefore,the research community has put little effort into the investigation of dynamics of systems with inertia variation.However,the mass change of the receiver is up to 18.8%within a period of 60 s in aerial refueling.Under this condition,the effect of the mass change on the dynamics of the aircraft could not be negligible.Earlier aerial refueling studies7,8either ignore the effect of mass transfer completely or treat it as a disturbance causing parametric uncertainty.To develop an accurate dynamic model of the aircraft,Dogan et al.9–11proposed a set of nonlinear,6-DOF equations of motion of a receiver aircraft undergoing aerial refueling.The equations,including the time-varying inertia associated with fuel transfer as well as the vortex-induced wind effect from the tanker,are expressed relative to the tanker and very rigorous from the mathematical perspective.Even so,they have a very complex mathematical expression and seem too complicated to be easily applied to the simulation and analysis.

    This paper focuses on the development of the derivation of a simpler aircraft model with the time-varying mass and inertia property.By extending equations of motion of a fixed weight aircraft,5comparing with the ones developed by Dogan et al.,the equations can reflect the time-varying mass and inertia property of an aircraft more factually due to(A)reflecting both the fuel increase and the fuel decrease at the same time(i.e.,the fuel transfer in refueling and the fuel consumption by engines),and(B)the fuel tanks having a rectangle shape rather than a mass point,(C)the equations of motion are derived relative to an inertial frame rather than the tanker body frame,(D)the model has a simpler mathematical expression and is easily applied to the simulation and analysis.Thus,they have a more strong universality.

    To validate the proposed aircraft model in an integrated simulation environment for probe-drogue-based autonomous aerial refueling,a dynamic model of the HDA developed by the authors12–16and a dynamic model of the vortex effect in aerial refueling developed by Dogan et al.17–19are also used in the simulation.

    2.Equations of motion of a time-varying inertia aircraft

    2.1.Equations of motion relative to Earth-centered inertial frame

    2.1.1.Definitions of coordinate frames and assumptions

    A receiver aircraft often performs long-range,even globestraddling military operations.Under this condition,the effect of the shape and the rotation of the Earth on the dynamics of the receiver aircraft should not be totally omitted.Therefore,without any loss of generality,an Earth-centered inertial(ECI)frame is chosen as an inertial reference frame.Fig.1 shows the ECI frame(Oexeyeze),an aircraft-body coordinate(ABC)frame(Obxbybzb),and an intermediate north-eastdown(NED)frame(Onxnynzn)on the surface of the Earth,5and ωeis the absolute angular velocity vector of Earth’rotation,rnis the position vector of the center of mass of the aircraft inOnxnynzn.

    To facilitate the derivation of the dynamics equations including the effect of time-varying inertia,an aircraft-fuel system is considered to comprise two parts as(A)solid,(B)fuel in tanks11as shown in Fig.1.The total amount of mass that occupies the fuel pipes has not been included here.The solid part is considered rigid with constant mass.The origin ofObxbybzbis chosen to be at the center of mass of the solid part,such that the standard definitions of aerodynamic variables and aerodynamic stability derivatives can be directly used without any modification or reinterpretation.5The fuel in fuel tanks is represented bykmasses.The mass of fuel in thejth fuel tank,mj(j=1,2,...,k),is time varying.The position vector from the center of mass of the fuel to the origin ofObxbybzb,rj,is assumed to be a constant.

    2.1.2.Translation kinematics

    Taking a derivative of the inertial position vector reof the center of mass of the receiver aircraft inOexeyezeyields5

    Fig.1 Coordinate frames for modeling.

    where ωe× redenotes the absolute velocity of surrounding air at re,B the rotation matrix fromOexeyezetoObxbybzb,BBthe rotation matrix fromOnxnynzntoObxbybzb,Vebthe absolute velocity vector of the aircraft center of mass expressed inObxbybzbrelative toOexeyeze,Vbthe relative velocity vector of the aircraft center of mass expressed inObxbybzbrelative toOexeyeze.φ,θ,ψ are the orientation angles.l,μ are the longitude and latitude.

    From Eq.(1),the translational kinematic equation of the aircraft relative toOexeyezecould be expressed as

    2.1.3.Translation dynamics

    Newton’s law for translational motion of a time-varying inertia aircraft takes the form

    where Fbis the vector sum of the aerodynamic and propulsion forces expressed inObxbybzb,mthe overall weights of the aircraft,gethe gravitational acceleration vector measured inOexeyeze.

    According to the theorem of Coriolis5,the differentiation of Vbcould be expressed as

    where ωbis the absolute angular velocity vector of the aircraft.

    It is straightforward to know˙ωe=0.Substituting Eqs.(5)and(7)into Eq.(6)yields

    wheremcould be regarded as the sum of the mass of the solid part and the fuel in all tanks and expressed as

    wheremsis the mass of the solid part of the receiver.

    It is straightforward to know˙ms=0.Thus,˙mcould be expressed as

    Substituting Eqs.(1),(9),and(10)into Eq.(8)and rearranging the equation,the translation dynamics equation of a receiver with varying mass relative toOexeyezecould be expressed as

    Known from Eq.(11),if the dynamics effects of varying inertia,i.e.,mjand˙mj,are neglected,Eq.(11)will have the same expression as a fixed weight aircraft’s translation dynamics equation.

    2.1.4.Rotational kinematics

    The rotational kinematical differential equations take the form of Poisson’s equation5

    where Ω is the cross-product matrix of the angular rates with respect to the axes of the rotating frame5,that is

    wherep,q,rare roll,pitch,and yaw angular rates respectively.

    2.1.5.Rotational dynamics

    The aerodynamic effects depend on the shape of the receiver aircraft body and the motion of the air relative to the body.Thus,the aerodynamic moments acting on a fixed weight aircraft can be directly used without any modification or reinterpretation here.To simplify the analysis,the thrust vector is aligned alongObxbybzb’sxaxis through the center of mass of the receiver,pointing toward its nose.However,the gravity vector of the fuel in each tank does not act through the center of mass of the receiver due to the tanks’configuration.Thus,the total moment T about the center of mass of the receiver could be expressed as

    Newton’s law for rotational motion takes the form

    where I is the total inertia matrix of the receiver.

    According to the theorem of Coriolis,5Eq.(15)could be written as

    where I could be regarded as the sum of the inertia matrix of the solid part and the fuel in all tanks,and expressed as follows

    where Isis the sum of the inertia matrix of the solid part of the receiver and Ijthe inertia matrix of the fuel in thejth fuel tank.

    According to the parallel-axis theorem,Ijcould be written as

    where I0jis the inertia matrix of the fuel in thejth fuel tank relative to its own parallel-axis and E the unit matrix.

    It is straightforward to know that˙Is=0.Thus,˙I could be expressed as

    Substituting Eqs.(14),(17)–(19)into Eq.(16),and rearranging the equation,the rotation dynamics equation of the receiver with varying inertia could be expressed as

    Known from Eq.(20),if the dynamics effects of varying inertia,i.e.,˙mjand˙I0j,are neglected,Eq.(20)will have the same expression as a fixed weight aircraft’s rotation dynamics equation.

    2.2.Equations of motion relative to NED frame

    For many short range flights,the effect of the shape and the rotation of the Earth on the dynamics of a receiver aircraft could be neglected.Therefore,equations of motion of a time-varying inertia receiver relative to the NED frame are more practical.Assuming ωe=0,replacing B with BB,and replacing gewith g,a reduced set of equations of motion of a time-varying inertia receiver relative to the NED frame is listed as follows:

    Navigation equation is

    Force equation is

    where g is the gravitational acceleration vector measured inOnxnynzn.

    Kinematic equation is

    Moment equation is

    Eq.(23),which consists of 9 coupled differential equations and is highly redundant,could be simplified by two methods.One is three-variable attitude(Euler angles)method commonly used for simulation,and the other is four-variable attitude method suitable for around-the-Earth flight,all-attitude flight and simulation of spinning bodies.

    2.3.Fuel tank configuration and fuel consumption

    A tailless fighter aircraft with innovative control effectors(ICE)11as the receiver is employed to show how the fuel tanks areincorporated asmathematicalquantities,which are required in the equations of motion.It is assumed that the amount and the location of the tanks,the fuel capacity of each tank,and the fuel flow distribution are set identical to Ref.11.However,each fuel tank has a rectangular shape rather than a lumped mass.Fig.2 shows the sketch of the ICE with 4 fuel tanks and thejth rectangular fuel tank with varying fuel.

    For a rectangular fuel tank as shown in Fig.2,I0jin Eq.(18)could be expressed as

    Fig.2 Configuration of fuel tanks of ICE.

    where the height of the fuelczcould be expressed as a function about the mass of the fuel in thejth tank as

    where ρ is the density of the fuel.ax,byare shown in Fig.2.

    For the rectangular fuel tank as shown in Fig.2,˙I could be expressed as

    The mass of the fuel in the fuel tanks will decrease due to the fuel consumption by engines.Thus,it is assumed that the consumption fuel flow by the engine˙mEcould be expressed as a function about the thrustT:

    whereamis the acceleration of the fuel consumption by the engine and η the equivalent coefficient.

    2.4.Position-tracking LQR controller for approaching and station-keeping in aerial refueling

    To ensure the safety and avoid the core area of the vortexinduced wind field,a reference trajectory is generated by a polynomial fitting curve20through the formation position(FP),the turning position(TP),the pre-docking position(PP),and the drogue position(DP)as shown in Fig.3.

    To employ the LQR control,the nonlinear equations of motion of the receiver aircraft relative to the NED frame derived above is trimmed and linearized at the steady level flight condition and initial mass.Linearizing Eqs.(21)-(24)yields:

    Fig.3 Position transition of ICE from approaching to docking.

    where x=x-x0,x=[V,α,β,p,q,r,ψ,θ,φ,xU,yU,zU]T,Vis the airspeed,α the angle of attack,β the angle of sideslip,xU,yU,zUare the positions;x0is the value of the trimmed states.u=u-u0,u=[δp,δe,δc,T],u0is the value of the trimmed inputs,δp,δe,δcare the angles of the surface.A0,B0are the state matrix and control matrix.

    To ensure zero tracking error at the steady state condition,three integrals of the position error as additional states are augmented into the state vector as

    The augmented state-space equation becomes

    Using LQR design technique,the state feedback gain matrix K is obtained to minimize the cost function:

    where Q is symmetric positive semi-definite and R symmetric positive definite.

    Selecting Q,R and solving the Riccatiequations PA+ATP-PBR-1BTP+Q=0 yields the input

    where K=R-1BTP,P is the symmetric positive semi-definite matrix.

    3.Simulation results

    3.1.Model of hose–drogue assembly with a nonlinear controller

    The HDA is modeled by a sequence of variable-length links connected with frictionless joints as illustrated in Fig.412–16,where all variables can be found in Refs.12–16A set of iterative equations of the hose’s three-dimensional motion is derived subject to hose reeling in/out,tanker motion,hose restoring force due to bending,gravity,and aerodynamic loads accounting for the effects of steady wind,atmospheric turbulence,tanker vortex.Two nonlinear control strategies are designed to suppress the hose whipping phenomenon.For details of the model and its controller,see the research previously conducted by the authors.12–16

    3.2.Modeling wind effects

    The aerodynamic forces and moments on a receiver aircraft are created by its motion relative to the surrounding air.5The air itself is in motion relative to the inertial frame to form the wind.The wind effects acting on the receiver aircraft should consist of two parts,the uniform local winds and the nonuniform wind induced by the wake vortex of the tanker.

    The effects of the uniform local winds,including the gust,the turbulence,the horizontal wind,and the wind shear,could be directly applied to the center of mass of the receiver aircraft.To highlight the mass change effects on the receiver’s states,verylightDryden wind turbulence(wind speed range:0–0.03 m/s)is carried out in simulation,because attitudes(φ,θ,ψ),angle of attack and side-slip angle(α,β)of the receiver in the close range aerial refueling scenario are little and serious Dryden wind turbulence will cover up them.

    However,the vortex-induced wind field acting on the receiver aircraft is nonuniform in nature,which could not be directly applied to the center of mass of the receiver aircraft.For simplicity,with a 39.88 m wing span and a 2.85 m distance between the pod and the top of the right wing,the vortex model and averaging technique derived by Dogan et al.17–19are employed to approximate the effects of the vortexinduced wind field on the center of mass of the receiver aircraft.A vertical slice of the wake 10 m behind the tanker is illustrated in Fig.5.The equivalent component of the wind vectorVwzacting on the receiver’s center of mass along the tanker’szaxis in the vortex-induced wind field is shown in Fig.6

    Thus,the velocity Vrand angular velocity vector ωrof the center of mass of the receiver aircraft with respect to the surrounding air are given by

    where Vwand ωware the velocity vector and angular velocity vector of the uniform local winds.Vvand ωvare the equivalent velocity vector and equivalent angular velocity vector induced by the tanker vortex.

    Then,to be able to use the equations of motion of the receiver aircraft without modification,Vrand ωr,rather than Vband ωb,could be used in the calculation of aerodynamic forces and moments.

    Fig.4 Configuration and modeling assumptions of a hose–drogue aerial refueling system.12–16

    Fig.5 Vertical slice of tanker wake 10 m behind tanker.

    3.3.Simulation of autonomous aerial refueling

    To validate the integrated simulation environment for probedrogue-based autonomous aerial refueling,the nonlinear equations of motion of the ICE relative toOnxnynznare implemented in simulation along with KC-135R.Both the receiver and the tanker models are equipped with their own controller based on LQR.The nominal altitude of the tanker is 7010 m with a speed of 200 m/s in a steady level flight.Other parameters used in simulation are set as follows.

    (1)The center of mass positions of the forward tanks and the aft tanks are(4,±4,0.2)and(-4,±4,0.2)respectively relative to the receiver aircraft body frame.a1x=a2x=b1y=b2y=b3y=b4y=2 m, and a3x=a4x=4 m.The initial quantity of the fuel in each tank is set to be 10%of the total capacity,i.e.,214.8 kg in Tank 1 and Tank 2,389.8 kg in Tank 3 and Tank 4.For simplicity,it is assumed that the engine pumps the identicalfuelflow distribution from each tank.am=0.12 kg/s2.Other parameters about ICE can be found in Ref.11

    (2)LOR controller gains: Q=10-2× diag{1,1,1,1,1,1,1,1,1,2,0.01,2,2,0.01,2},R=10-3×diag{0.00001,1,1,1}.

    (3)Relative to the tanker’s body frame,FP,TP,and PP are set equal to(-31,34,2),(-52.8,34,12),and(-52.8,18.2,12).DP is a variable depending on the position of the drogue,which is decided by the HDA model in previous research conducted by the authors.12-16

    Fig.6 Equivalent component of wind vector along tanker’s z axis in the vortex-induced wind field.

    3.3.1.Approaching and docking phase

    The ICE is initially at(-31,34,2)m inObxbybzbof the tanker and should approach and dock with a drogue for refueling along the reference trajectory programmed above.Figs.7–13 depict the results of the simulations when executed for an approaching and docking flight before refueling.

    Fig.7 shows the actual trajectory and the reference trajectory of the ICE relative toObxbybzbof the tanker during approaching and docking.Fig.8 illustrates the triaxial details and the tracking errors of the trajectory.As seen from Figs.7 and 8,under the influences of the fuel consumption and the aerodynamic disturbance,the LQR controller still enables the ICE to track the reference trajectory exactly.

    Figs.9 and 10 show the time histories of the mass of the fuel and the inertia of the ICE due to fuel consumption by the engine during approaching and docking.

    As seen from Fig.9 and Fig.10,the mass of the fuel and the inertia of the ICE decrease due to fuel consumption by the engine throughout the approaching and docking flight.As a result of the fuel tank symmetry,Ixy,Iyx,Iyz,andIzyalways maintain to be 0.From Eq.(28),the rate of the fuel consumption by the engine depends on the thrust.Therefore,the rate of the fuel consumption increases gradually in pace with the transition from the deceleration to the acceleration of the receiver as shown in Fig.11.

    The true airspeed,the aerodynamic angles,and the orientation of the ICE are presented in Figs.11 and 12.Fig.13 illustrates the values of the elevon δe,pitch flap δp,and clamshell δcdeflections and the thrust.

    As shown in Figs.5–7,outside the wing of the tanker,the local airflow of the vortex flows upward.If the controller does not respond,the angle of attack α will increase,thexcomponent of the true airspeedVwill decrease,and then the receiver will not be able to track the reference trajectory exactly.Therefore,from the FP to the TP,about 0–150 s,the controller negatively increases δpto decrease the pitch angle θ of the ICE to resist continuous growth of α.Thexcomponent of the gravity of the ICE due to the decreased pitch angle,except for the effect of the deceleration of the ICE,will decrease the thrust.As a result,the fuel consumption will be decreased with respect to a level flight.However,the condition inside the wing of the tanker,about 150–300 s,is exactly opposite of a bond.Due to the fuel tank symmetry,lateral states and control variables in the steady state tend to 0 in the presence of the disturbance caused by the vortex and the Dryden wind turbulence.

    Fig.7 Actual trajectory and reference trajectory of ICE relative to tanker.

    Fig.8 Triaxial details and tracking errors of trajectory.

    Fig.9 Time histories of fuel mass considering fuel consumption by engine.

    Fig.10 Time histories of inertia of ICE considering fuel consumption by engine.

    Fig.11 Time histories of true airspeed and aerodynamic angles.

    Fig.12 Time histories of orientation.

    3.3.2.Station-keeping and refueling phase

    After docking,the ICE should maintain a steady position relative toObxbybzbof the tanker to receive 0.04416 m3/s of JP-4 fuel.Figs.14–19 depict the results of the simulations during the station-keeping and refueling phase for three different cases as the same as three cases used in Ref.11Case 1:both forward tanks refuel first,but at half of the maximum flow rate,until each is half full.Next,both aft tanks refuel,again at half of the maximum flow rate,until each is half full.Case 2:neither of the right-hand fuel tanks receives fuel,but the lefthand forward tank is fully fueled first,followed by the lefthand aft.Case 3:the right-hand forward and left-hand aft tanks do not fill,but the other tanks fully fuel,with the forward tank refueling first.To allow for better comparisons,the fuel tanks in Case 1 only fill to half capacity so that the transition from refueling the forward tanks to the aft tanks occurs at the same time in all three cases.In addition,the postrefuel total aircraft weight is identical for three cases.In each case,refueling starts att=300 s.

    Figs.14 and 15 illustrate the ICE’s position with respect to the docking position.The controller,as expected,has the same station-keeping performance for three cases except spikes inxandzposition deviations when fuel flow starts.The variations of longitudinal states of the ICE,i.e.,thexandzdeviations from the docking position,V,α,θ,δp,andTas seen later,have the same tendencies.

    Figs.9 and 10 show the time histories of the mass of the fuel and the inertia of the ICE due to the fuel transfer as well as the fuel consumption by engines during station-keeping and refueling.Each fuel tank is filled with the fuel as stated above for three cases.To the end of refueling at about 467 s,the total fuel filled in the ICE(eliminating the fuel consumed by engines)is 5670.5 kg accounting for 46.06%of the initial weight of the ICE.

    As shown in Figs.16 and 17,Ix,Iy,Iz,Ixz,andIzxhave identical tendencies for all three cases due to thexOzplane of symmetry of the ICE.Under this condition,for every productyjzjoryjxjin an inertia computation there is a different value for a different case.Therefore,Ixy,Iyx,Iyz,andIzyare different for three cases.

    Fig.13 Time histories of actuators and thrust.

    Fig.14 ICE’s position variations with respect to docking position.

    Fig.15 Triaxial details of ICE’s position variations with respect to docking position.

    Fig.16 Time histories of fuel mass during fuel transfer.

    Fig.17 Time histories of inertia of the ICE during fuel transfer.

    Fig.18 Time histories of orientation of the ICE during fuel transfer.

    Fig.19 Time histories of actuators and thrust of ICE during fuel transfer.

    Fig.18 illustrates the deviations of the airspeed and the orientation respectively.The controller,as expected,has the same speed-keeping performance for three cases except spikes when fuel flow starts.As stated above,the variations of longitudinal states of the ICE for three cases have the same tendencies.Before aft tanks are filled with the fuel for three cases,α and θ are increased by the controller to increase the upward moment to compensate for the effects of the forward tank.At this time,the deviations of the lateral states are still small.Therefore,α and θ have identical trajectories for all three cases.After the aft tanks are filled with the fuel for three cases,the gravity of the aft tanks grows gradually and lowers the values of the longitudinal states than before.For Cases 2 and 3,the effects caused by the fuel transfer could be regarded as identical for the longitudinal states of the ICE based on thexOzplane of symmetry.However,the asymmetric refueling in Cases 2 and 3 produces an asymmetric moment in the lateral direction due to the fuel tank gravity.To keep lateral balance and track the commanded speed,the controller makes the ICE roll or yaw a small angle in the corresponding direction to compensate for the opposite-side tank refueling.As a result,the component of the lift in the gravitational direction decreases and must be compensated by larger longitudinal states than Case 1.When compared to the same side forward to aft tank transition of Case 2,the opposite-direction tendencies of the sideslip,roll,and yaw angles immediately follow the Case 3 transition from the forward to the diagonal aft fuel tank to compensate for the opposite-side tank refueling and maintain a commanded steady level flight.After completion of the refueling,the longitudinal trim states are a little larger than their initial values to produce enough lift.The roll and yaw angles are the smallest for Case 1,which results in zero deviations as expected based on the symmetry,while the roll and yaw trim angles are different for Cases 2 and 3.This demonstrates the effect of the asymmetric mass distribution on the steady-state trim condition.

    Fig.19 illustrates the values of δe,δp,δcdeflections and the thrust for three cases.A comparison can be made between Figs.18 and 19 to show how the effectors react to the refueling in three cases and how their steady-state positions correspond to the receiver’s steady-state Euler angles.The longitudinal control variables,i.e.,δpand the thrust,match the longitudinal states shown in Fig.18.When the refueling of Case 3 transitions from the forward to the diagonal aft tank,both δeand δcreact immediately to compensate for the opposite-side tank refueling.The magnitude of the effecter deflections changes steadily throughout the refueling to account for the steadily growing fuel tank masses.This is caused by the contribution of the gravitational moment caused by the fuel masses to the total moment of the external forces about the origin ofObxbybzbof the ICE.The thrusts in Fig.19 for three cases,which prove to be very similar,tend to grow as the refueling occurs to keep the receiver at the commanded velocity while experiencing mass increase.The steady-state thrust value is somewhat higher than initially due to the higher postrefuel aircraft weight.

    For ease of comparison,two of the simulation results in Ref.11are cited as Figs.20 and 21.It is known to all that if the fuel mass in one tank is refueled gradually,the states of the receiver should also change gradually with it.The tendency of the state variation should not change unless the refueling course changes.But known from Figs.20 and 21,there exists an obvious state mutation(especially the yaw and the clamshell),when each thank is filled,which goes against the continuity of the physical system as stated above.Taking the clamshell in the Case 3 as an example,the deflection should monotonically increase to the maximum in a direction.Then the deflection should go through the same process but in the other direction,when the transition from refueling the forward tanks to the aft tanks occurs,because the fuel is always filled at a constant flow rate.Besides,the deflection should always maintain a maximum in the relevant direction after the refueling stops.Therefore,although the model in Ref.11have a correct mathematical derivation,some mistakes occur in the simulation application.Fortunately,this phenomenon does not appear at all in the model proposed here.

    Fig.20 Time history of receiver relative orientation deviation in straight level flight.

    Fig.21 Time history of receiver control surface deflections in straight level flight.

    4.Conclusions

    (1)Taking into account three major factors,i.e.,the dynamic models of the receiver with the time-varying inertia property,an HDA dynamic model with the variable-length property,and the wind effect due to the tanker’s trailing vortices,an integrated simulation environment for probe-drogue-based autonomous aerial refueling is developed.The environment could be used to further develop and evaluate aircraft control,sensor systems,dynamics and modeling in the autonomous aerial refueling research.

    (2)Considering a rectangular shape of the tanks rather than a lumped mass,a set of new simpler equations of motion of an aircraft with the time-varying inertia property due to fuel transfer as well as the fuel consumption by engines is derived relative to the ECI frame or the NED frame.

    (3)The ICE as the receiver and a LQR-based autonomous flight controller is employed to validate the integrated simulation environment for probe-drogue-based autonomous aerial refueling.

    Acknowledgments

    This research was supported by the National Natural Science Foundation of China(Nos.61473307,61304120).

    1.Thomas PR,Bhandari U,Bullock S,Richardson TS.Advances in air to air refueling.Prog Aerosp Sci2014;71:14–35.

    2.Quan Q,Wei ZB,Gao J,Zhang RF,Cai KY.A survey on modeling and control problems for probe&drogue autonomous aerial docking.Acta Aeronaut Astronaut Sin2014;35(9):2390–410[Chinese].

    3.Lu YP,Yang CX,Liu YY.A survey of modeling and control technologies for aerial refueling system.Acta Aeronaut Astronaut Sin2014;35(9):2375–89[Chinese].

    4.Dong XM,Xu YJ,Chen B.Progress and challenges in automatic aerial refueling.J Air Force Eng Univ(Nat Sci Ed)2008;9(6):1–5[Chinese].

    5.Lewis BL,Stevens FL.Aircraft control and simulation.New York:Wiley;1992.p.21.

    6.Roskam J.Airplane flight dynamics and automatic flight controls.Lawrence:DAR Corporation;2001.p.5.

    7.Bennington MA,Visser KD.Aerial refueling implications for commercial aviation.J Aircraft2005;42(2):366–75.

    8.Bloy A,Khan M.Modeling of the receiver aircraft in air-to-air refueling.J Aircraft2001;38(2):393–6.

    9.Venkataramanan S,Dogan A.Dynamic effects of trailing vortex with turbulence&time-varying inertia in aerial refueling.AIAA atmospheric flight mechanics conference.Reston:AIAA;2004.

    10.Tucker J,Dogan A.Derivation of the dynamics equations of receiver aircraft in aerial refueling.45th AIAA aerospace sciences meeting.Reston:AIAA;2007.

    11.Waishek J,Dogan A,Blake W.Derivation of the dynamics equations of receiver aircraft in aerial refueling.J Guid Control Dyn2009;32(2):585–97.

    12.Wang HT,Dong XM,Dou HF,Xue JP.Dynamic modeling and characteristics analysis of hose-paradrogue aerial refueling systems.J Beijing Univ Aeronaut Astronaut2014;40(1):92–8[Chinese].

    13.Wang HT,Dong XM,Xue JP,Liu JL.Dynamic modeling of hose–drogue aerial refueling system and integral sliding mode backstepping control for whipping phenomenon.Chin J Aeronaut2014;27(4):930–46.

    14.Wang HT,Dong XM,Liu JL,Wang J.Dynamics and control of the hose whipping phenomenon in aerial refueling.2015 IEEE aerospace conference.Piscataway,NJ:IEEE Press;2015.

    15.Wang HT,Dong XM,Guo J,Liu JL,Wang J.Modeling and analysis of the hose whipping phenomenon of the refueling hose–drogue assembly.Acta Aeronaut Astronaut Sin2015;36(9):3116–27[Chinese].

    16.Wang HT,Dong XM,Guo J,Liu JL,Wang J.Modeling and simulation for autonomous aerial refueling using a probe-drogue system.AIAAatmosphericflightmechanicsconference.Reston:AIAA;2015.

    17.Dogan A,Kim E,Blake W.Control and simulation of relative motion for aerial refueling in racetrack maneuvers.J Guid Control Dyn2007;35(5):1551–7.

    18.Dogan A,Venkataramanan S.Nonlinear control for reconfiguration of unmanned-aerial-vehicle formation.J Guid Control Dyn2005;28(4):667–78.

    19.Dogan A,Venkataramanan S,Blake W.Modeling of aerodynamic coupling between aircraft in close proximity.J Aircraft2005;42(4):941–55.

    20.Tandale MD,Bowers R,Valasek J.Trajectory tracking controller for vision-based probe and drogue autonomous aerial refueling.J Guid Control Dyn2006;29(4):846–57.

    21 January 2015;revised 19 May 2015;accepted 11 August 2015

    Available online 5 March 2016

    ?2016 Chinese Society of Aeronautics and Astronautics.Published by Elsevier Ltd.This is an open access article under the CC BY-NC-ND license(http://creativecommons.org/licenses/by-nc-nd/4.0/).

    *Corresponding author.Tel.:+86 29 84787400 801.

    E-mail addresses: wanghaitao198638@163.com (H. Wang),dongxinmin@139.com(X.Dong).

    Peer review under responsibility of Editorial Committee of CJA.

    Wang Haitaoreceived his B.S.and M.S.degrees from Air Force Engineering University in 2009 and 2011,respectively.He is now pursing his Ph.D.degree there.His main research interests include system modeling,control and simulation of aerial refueling.

    Dong Xinminis a professor at Air Force Engineering University.He received his B.S.and M.S.degrees from Northwestern Polytechnical University in 1983 and 1988,respectively,and then received his Ph.D.degree from Xi’an Jiaotong University in 1991.He has worked in the field of flight control for more 20 years,and his current research interests include modern control theory,machine vision,and autonomous aerial refueling.

    Xue Jianpingis an associate professor at Air Force Engineering University.He received his B.S.and M.S.degrees from Air Force Engineering University in 1992 and 2005,respectively.His current research interests include computer control and system simulation.

    Liu Jiaolongreceived his B.S.and M.S.degrees from Air Force Engineering University in 2011 and 2013,respectively.He is now pursing his Ph.D.degree there.His main research interests include CFD analysis and learning control.

    Wang Jianreceived his B.S.degree from Air Force Engineering University in 2013.He is now pursing his M.S.degree there.His main research interests include CFD analysis,control and simulation of aerial refueling.

    欧美另类亚洲清纯唯美| 亚洲中文字幕日韩| 精品久久久久久久人妻蜜臀av | 亚洲第一欧美日韩一区二区三区| av超薄肉色丝袜交足视频| 久久中文字幕人妻熟女| 国产精品久久久久久精品电影 | 午夜视频精品福利| 少妇 在线观看| 久久青草综合色| 97超级碰碰碰精品色视频在线观看| 亚洲aⅴ乱码一区二区在线播放 | 久久这里只有精品19| 亚洲aⅴ乱码一区二区在线播放 | 一边摸一边做爽爽视频免费| 精品午夜福利视频在线观看一区| 好男人在线观看高清免费视频 | 国产av精品麻豆| 99国产综合亚洲精品| 涩涩av久久男人的天堂| 欧美成人午夜精品| 久久青草综合色| 亚洲av电影不卡..在线观看| 这个男人来自地球电影免费观看| 变态另类丝袜制服| 日韩 欧美 亚洲 中文字幕| 日韩欧美国产一区二区入口| 亚洲中文字幕一区二区三区有码在线看 | 精品人妻1区二区| 免费在线观看影片大全网站| 亚洲欧美精品综合一区二区三区| 丝袜美足系列| 99国产精品99久久久久| 久久精品亚洲熟妇少妇任你| 黄网站色视频无遮挡免费观看| 午夜精品国产一区二区电影| 女人高潮潮喷娇喘18禁视频| 国产精品久久久久久亚洲av鲁大| 久久久国产欧美日韩av| 亚洲国产毛片av蜜桃av| 久久天堂一区二区三区四区| 亚洲激情在线av| 亚洲视频免费观看视频| 久久狼人影院| 波多野结衣高清无吗| 久久久国产欧美日韩av| 午夜视频精品福利| 久久久久久人人人人人| 久久精品91蜜桃| 侵犯人妻中文字幕一二三四区| 免费观看精品视频网站| 日本vs欧美在线观看视频| 亚洲色图 男人天堂 中文字幕| 日日干狠狠操夜夜爽| 久久中文字幕人妻熟女| 国产亚洲欧美在线一区二区| 欧美黑人精品巨大| 人人妻,人人澡人人爽秒播| 麻豆成人av在线观看| 在线国产一区二区在线| 可以在线观看毛片的网站| 国产高清videossex| 桃色一区二区三区在线观看| 18禁黄网站禁片午夜丰满| 国产高清有码在线观看视频 | 男女之事视频高清在线观看| 国产亚洲欧美在线一区二区| 久久精品人人爽人人爽视色| 亚洲激情在线av| 午夜免费成人在线视频| 在线观看www视频免费| 欧美精品亚洲一区二区| 亚洲国产日韩欧美精品在线观看 | 91麻豆精品激情在线观看国产| 国产男靠女视频免费网站| 色哟哟哟哟哟哟| 国语自产精品视频在线第100页| 丝袜美足系列| 国产激情欧美一区二区| 精品国产国语对白av| 淫妇啪啪啪对白视频| 午夜成年电影在线免费观看| 国产成人av教育| 欧美日韩一级在线毛片| 中出人妻视频一区二区| 亚洲avbb在线观看| 精品欧美一区二区三区在线| 亚洲欧美日韩高清在线视频| 多毛熟女@视频| 在线播放国产精品三级| 国产99久久九九免费精品| 成人精品一区二区免费| 色婷婷久久久亚洲欧美| 亚洲中文字幕一区二区三区有码在线看 | 免费在线观看黄色视频的| 桃红色精品国产亚洲av| 国产成人av教育| 精品日产1卡2卡| www.精华液| 成人亚洲精品av一区二区| 级片在线观看| 久久伊人香网站| 国内精品久久久久精免费| 国产人伦9x9x在线观看| 久久国产精品人妻蜜桃| av天堂久久9| 啦啦啦韩国在线观看视频| 女性被躁到高潮视频| 亚洲 欧美 日韩 在线 免费| 黄色 视频免费看| 成人特级黄色片久久久久久久| 国产亚洲精品久久久久久毛片| 亚洲熟女毛片儿| av天堂在线播放| 国产真人三级小视频在线观看| 国产主播在线观看一区二区| 久久精品aⅴ一区二区三区四区| 波多野结衣av一区二区av| 啦啦啦免费观看视频1| 一区福利在线观看| 国产精品香港三级国产av潘金莲| 久久精品aⅴ一区二区三区四区| 美女高潮喷水抽搐中文字幕| 欧美日本亚洲视频在线播放| av欧美777| 天天躁狠狠躁夜夜躁狠狠躁| 欧美日本亚洲视频在线播放| 好男人电影高清在线观看| 亚洲精品国产区一区二| av免费在线观看网站| 国产午夜福利久久久久久| 日韩免费av在线播放| 午夜福利高清视频| 成人av一区二区三区在线看| 欧美一区二区精品小视频在线| avwww免费| 宅男免费午夜| 视频区欧美日本亚洲| 精品免费久久久久久久清纯| 久久精品国产亚洲av高清一级| www.www免费av| 国产极品粉嫩免费观看在线| 国内精品久久久久久久电影| 99久久精品国产亚洲精品| 欧美日韩亚洲国产一区二区在线观看| 国产av一区在线观看免费| 中文字幕最新亚洲高清| 日本在线视频免费播放| 成人18禁在线播放| 国内久久婷婷六月综合欲色啪| 亚洲第一电影网av| 真人做人爱边吃奶动态| 国产视频一区二区在线看| 亚洲伊人色综图| 欧美黑人欧美精品刺激| 欧美日韩亚洲综合一区二区三区_| 久久性视频一级片| 性欧美人与动物交配| 欧美大码av| 成人永久免费在线观看视频| 不卡一级毛片| 巨乳人妻的诱惑在线观看| 国产欧美日韩一区二区三| 国产精品二区激情视频| 黄色视频,在线免费观看| 国产成人影院久久av| 国产精品秋霞免费鲁丝片| 精品日产1卡2卡| 亚洲中文日韩欧美视频| 亚洲一区二区三区不卡视频| 欧美精品亚洲一区二区| 午夜视频精品福利| 国产精品久久久久久精品电影 | 黑人巨大精品欧美一区二区蜜桃| 在线视频色国产色| 自线自在国产av| 纯流量卡能插随身wifi吗| 嫩草影视91久久| 欧美亚洲日本最大视频资源| 久久欧美精品欧美久久欧美| 在线视频色国产色| 亚洲精品久久成人aⅴ小说| 岛国视频午夜一区免费看| 亚洲性夜色夜夜综合| 日韩av在线大香蕉| 好男人在线观看高清免费视频 | 国产精品免费一区二区三区在线| 欧美黑人精品巨大| 亚洲自拍偷在线| 欧美亚洲日本最大视频资源| 国产精品亚洲美女久久久| 国产av一区在线观看免费| 久久青草综合色| 少妇熟女aⅴ在线视频| 国产精品亚洲一级av第二区| 国产成人欧美| 国产av在哪里看| 亚洲国产看品久久| 91精品三级在线观看| 免费高清在线观看日韩| 性色av乱码一区二区三区2| 日韩av在线大香蕉| 亚洲免费av在线视频| 国产精品一区二区在线不卡| 一本综合久久免费| 国产视频一区二区在线看| 欧美一区二区精品小视频在线| 欧美日本亚洲视频在线播放| 午夜福利一区二区在线看| 国产亚洲欧美在线一区二区| 久久久久久久精品吃奶| 女人被狂操c到高潮| 国产精品一区二区免费欧美| 国产亚洲精品一区二区www| 黄网站色视频无遮挡免费观看| 欧美日本亚洲视频在线播放| 欧美黑人欧美精品刺激| 国产伦人伦偷精品视频| 又大又爽又粗| aaaaa片日本免费| 亚洲天堂国产精品一区在线| 亚洲欧洲精品一区二区精品久久久| 韩国av一区二区三区四区| 国产视频一区二区在线看| 久久久久亚洲av毛片大全| 亚洲视频免费观看视频| videosex国产| 人成视频在线观看免费观看| 制服诱惑二区| 中文字幕色久视频| 黄色a级毛片大全视频| 午夜日韩欧美国产| 一级黄色大片毛片| 好男人电影高清在线观看| 久久国产精品男人的天堂亚洲| 国产xxxxx性猛交| 亚洲久久久国产精品| 18禁美女被吸乳视频| 香蕉国产在线看| 国产激情久久老熟女| 精品卡一卡二卡四卡免费| 不卡一级毛片| 亚洲视频免费观看视频| 免费看美女性在线毛片视频| 午夜精品在线福利| 亚洲av电影不卡..在线观看| 国产精品野战在线观看| 午夜两性在线视频| 亚洲精品在线观看二区| 岛国视频午夜一区免费看| 国产伦人伦偷精品视频| 国产精品亚洲av一区麻豆| 在线国产一区二区在线| 午夜精品国产一区二区电影| www.999成人在线观看| a在线观看视频网站| 熟女少妇亚洲综合色aaa.| 一本综合久久免费| 精品久久久久久,| 90打野战视频偷拍视频| 精品无人区乱码1区二区| 免费一级毛片在线播放高清视频 | 美女高潮到喷水免费观看| 欧美在线一区亚洲| 狂野欧美激情性xxxx| 一本大道久久a久久精品| 亚洲在线自拍视频| 少妇熟女aⅴ在线视频| 亚洲专区中文字幕在线| 他把我摸到了高潮在线观看| av视频在线观看入口| 亚洲avbb在线观看| 少妇熟女aⅴ在线视频| 亚洲国产精品sss在线观看| 久久人妻福利社区极品人妻图片| 亚洲成av人片免费观看| 精品久久久精品久久久| 久久久久国产精品人妻aⅴ院| 最好的美女福利视频网| 免费观看人在逋| 国产精品美女特级片免费视频播放器 | 此物有八面人人有两片| 欧美日本中文国产一区发布| 99香蕉大伊视频| 国产成人免费无遮挡视频| 精品不卡国产一区二区三区| 国产又色又爽无遮挡免费看| 满18在线观看网站| 国产麻豆69| 精品第一国产精品| 精品福利观看| 国产精品香港三级国产av潘金莲| 18禁美女被吸乳视频| 国产97色在线日韩免费| 在线观看舔阴道视频| 身体一侧抽搐| 亚洲精品国产精品久久久不卡| 国产色视频综合| 黄色a级毛片大全视频| 亚洲国产高清在线一区二区三 | 丝袜美足系列| 亚洲精品在线美女| 长腿黑丝高跟| 欧美日韩福利视频一区二区| 久久国产亚洲av麻豆专区| 黄色 视频免费看| 丰满的人妻完整版| 好男人在线观看高清免费视频 | 午夜久久久久精精品| 欧美精品亚洲一区二区| 欧美日本中文国产一区发布| 丝袜人妻中文字幕| 欧美中文综合在线视频| 午夜激情av网站| 国产又爽黄色视频| 亚洲七黄色美女视频| 可以免费在线观看a视频的电影网站| 看黄色毛片网站| 男女午夜视频在线观看| 亚洲狠狠婷婷综合久久图片| 嫁个100分男人电影在线观看| 成人18禁在线播放| 又紧又爽又黄一区二区| 国产欧美日韩一区二区三| 欧美成人一区二区免费高清观看 | 在线播放国产精品三级| 精品国产亚洲在线| 老汉色av国产亚洲站长工具| www.www免费av| 国产1区2区3区精品| 每晚都被弄得嗷嗷叫到高潮| 国产黄a三级三级三级人| 国产欧美日韩一区二区精品| 精品国内亚洲2022精品成人| 大型黄色视频在线免费观看| 亚洲av电影在线进入| 两人在一起打扑克的视频| 国产精品野战在线观看| 久久久久久人人人人人| 国产欧美日韩精品亚洲av| 一本久久中文字幕| 91老司机精品| 日韩中文字幕欧美一区二区| av天堂在线播放| 亚洲九九香蕉| 一进一出抽搐动态| 亚洲一区中文字幕在线| 欧美黑人精品巨大| bbb黄色大片| 日韩精品中文字幕看吧| 午夜福利,免费看| 操出白浆在线播放| 中文字幕高清在线视频| 午夜精品国产一区二区电影| 久久精品人人爽人人爽视色| 国内精品久久久久久久电影| 精品高清国产在线一区| 亚洲伊人色综图| 一级片免费观看大全| 亚洲欧洲精品一区二区精品久久久| 精品不卡国产一区二区三区| 亚洲久久久国产精品| 亚洲男人的天堂狠狠| 久久久久亚洲av毛片大全| 久久久久九九精品影院| 又黄又粗又硬又大视频| 免费在线观看完整版高清| 69精品国产乱码久久久| 一级毛片精品| 欧美日韩瑟瑟在线播放| 嫩草影视91久久| 制服丝袜大香蕉在线| 一卡2卡三卡四卡精品乱码亚洲| 精品国产亚洲在线| 91精品三级在线观看| 午夜福利在线观看吧| 免费在线观看日本一区| 一区二区三区高清视频在线| 亚洲片人在线观看| 18美女黄网站色大片免费观看| 久久人妻av系列| 精品人妻1区二区| 免费观看人在逋| 多毛熟女@视频| 妹子高潮喷水视频| 日韩三级视频一区二区三区| 黄色片一级片一级黄色片| 黑人操中国人逼视频| 精品久久久久久,| 午夜福利视频1000在线观看 | 久久精品aⅴ一区二区三区四区| 69精品国产乱码久久久| 亚洲专区字幕在线| 村上凉子中文字幕在线| 国产精品久久久久久人妻精品电影| 正在播放国产对白刺激| 伊人久久大香线蕉亚洲五| 亚洲性夜色夜夜综合| 多毛熟女@视频| 日韩精品青青久久久久久| 中文字幕另类日韩欧美亚洲嫩草| 亚洲人成电影免费在线| aaaaa片日本免费| 欧美乱码精品一区二区三区| 一区福利在线观看| 桃红色精品国产亚洲av| 色播在线永久视频| 国产精品久久久av美女十八| av超薄肉色丝袜交足视频| 窝窝影院91人妻| 日韩精品免费视频一区二区三区| 香蕉国产在线看| 咕卡用的链子| 国产高清激情床上av| 一本久久中文字幕| 一级毛片高清免费大全| 亚洲欧美日韩无卡精品| 最近最新免费中文字幕在线| 久久国产精品人妻蜜桃| 精品电影一区二区在线| 欧美日韩乱码在线| 日本三级黄在线观看| 亚洲成人免费电影在线观看| 国产欧美日韩精品亚洲av| 久久久精品国产亚洲av高清涩受| 亚洲自偷自拍图片 自拍| 丰满的人妻完整版| 久9热在线精品视频| 美女国产高潮福利片在线看| 亚洲,欧美精品.| 亚洲免费av在线视频| 在线观看免费视频日本深夜| 国产一区二区激情短视频| 自拍欧美九色日韩亚洲蝌蚪91| 国产精品免费视频内射| 久久香蕉国产精品| 最新在线观看一区二区三区| 久久国产亚洲av麻豆专区| 好男人电影高清在线观看| 精品免费久久久久久久清纯| 怎么达到女性高潮| 亚洲成人精品中文字幕电影| 国产一区二区三区综合在线观看| 精品无人区乱码1区二区| 欧美中文综合在线视频| av电影中文网址| 亚洲,欧美精品.| 女人被狂操c到高潮| 免费高清视频大片| 嫩草影院精品99| 欧美成人免费av一区二区三区| 亚洲欧美激情在线| 91成年电影在线观看| 校园春色视频在线观看| 亚洲欧美激情在线| 国产一区二区三区在线臀色熟女| 色综合欧美亚洲国产小说| 午夜两性在线视频| 91成年电影在线观看| 免费久久久久久久精品成人欧美视频| 国产成人欧美| a级毛片在线看网站| 亚洲av成人一区二区三| 两个人免费观看高清视频| 操出白浆在线播放| 欧美性长视频在线观看| 天堂动漫精品| 国产单亲对白刺激| 久久精品91蜜桃| 俄罗斯特黄特色一大片| 亚洲欧洲精品一区二区精品久久久| 一级a爱片免费观看的视频| 精品熟女少妇八av免费久了| 亚洲成人免费电影在线观看| 母亲3免费完整高清在线观看| 欧美日本视频| 欧美老熟妇乱子伦牲交| 久久草成人影院| 欧美一区二区精品小视频在线| 18禁裸乳无遮挡免费网站照片 | 欧美激情高清一区二区三区| 黄色a级毛片大全视频| av在线播放免费不卡| 在线观看免费午夜福利视频| 亚洲精品一卡2卡三卡4卡5卡| 亚洲av电影在线进入| 亚洲天堂国产精品一区在线| 午夜福利18| 国产国语露脸激情在线看| 不卡一级毛片| 大型av网站在线播放| 中文字幕高清在线视频| 亚洲第一青青草原| 桃色一区二区三区在线观看| 亚洲中文字幕一区二区三区有码在线看 | 精品久久久久久久人妻蜜臀av | 国产成人精品在线电影| 亚洲第一电影网av| 亚洲精品一区av在线观看| 久久精品91无色码中文字幕| 亚洲av片天天在线观看| 91麻豆精品激情在线观看国产| 国产亚洲欧美98| 纯流量卡能插随身wifi吗| 十八禁人妻一区二区| 欧美不卡视频在线免费观看 | 日韩精品中文字幕看吧| 国产精品自产拍在线观看55亚洲| 看免费av毛片| 男女床上黄色一级片免费看| 久久久久精品国产欧美久久久| 波多野结衣av一区二区av| 免费在线观看黄色视频的| 夜夜看夜夜爽夜夜摸| 亚洲精品久久国产高清桃花| 欧美国产精品va在线观看不卡| 90打野战视频偷拍视频| 真人做人爱边吃奶动态| 88av欧美| av天堂久久9| 精品国产美女av久久久久小说| 在线播放国产精品三级| 男人舔女人下体高潮全视频| 亚洲avbb在线观看| av欧美777| 色综合欧美亚洲国产小说| 久久国产精品男人的天堂亚洲| 三级毛片av免费| 亚洲专区中文字幕在线| 如日韩欧美国产精品一区二区三区| 欧美成人午夜精品| 日韩av在线大香蕉| 老司机深夜福利视频在线观看| 国产av精品麻豆| 精品久久久精品久久久| 欧美日韩中文字幕国产精品一区二区三区 | 欧美乱妇无乱码| 少妇被粗大的猛进出69影院| 国产精品秋霞免费鲁丝片| 日本免费一区二区三区高清不卡 | 咕卡用的链子| 日韩欧美免费精品| 日韩国内少妇激情av| 国产成人精品无人区| 亚洲五月色婷婷综合| 嫁个100分男人电影在线观看| 999久久久精品免费观看国产| 香蕉国产在线看| 久久婷婷成人综合色麻豆| 这个男人来自地球电影免费观看| 欧美精品亚洲一区二区| 亚洲av第一区精品v没综合| av电影中文网址| 天天添夜夜摸| 麻豆一二三区av精品| 长腿黑丝高跟| 欧美最黄视频在线播放免费| 国产精品九九99| 啦啦啦韩国在线观看视频| 亚洲第一青青草原| 国产精品永久免费网站| 午夜a级毛片| 国产片内射在线| 久久婷婷人人爽人人干人人爱 | 一区二区三区高清视频在线| 制服诱惑二区| 亚洲中文字幕一区二区三区有码在线看 | 亚洲自拍偷在线| 久久久久久亚洲精品国产蜜桃av| 亚洲专区中文字幕在线| 99精品久久久久人妻精品| 在线观看www视频免费| 免费看十八禁软件| 亚洲欧美激情在线| 午夜精品在线福利| 不卡av一区二区三区| 亚洲av熟女| 在线观看66精品国产| 日韩大码丰满熟妇| 国产视频一区二区在线看| a在线观看视频网站| 精品一区二区三区视频在线观看免费| 99riav亚洲国产免费| 国产免费男女视频| 中文字幕另类日韩欧美亚洲嫩草| 亚洲av片天天在线观看| 亚洲人成77777在线视频| 99国产精品99久久久久| 女性生殖器流出的白浆| 看片在线看免费视频| 淫秽高清视频在线观看| 午夜福利影视在线免费观看| 叶爱在线成人免费视频播放| 91九色精品人成在线观看| 成在线人永久免费视频| 搡老熟女国产l中国老女人| 久久久久久大精品| 热re99久久国产66热| 国产一区二区三区在线臀色熟女| 国产熟女xx| 久久中文字幕一级| 日韩国内少妇激情av| 国产成人精品无人区| 一边摸一边抽搐一进一出视频| 亚洲av熟女| 欧美成人午夜精品| 久久午夜亚洲精品久久| 午夜视频精品福利| 免费在线观看亚洲国产| 人人妻,人人澡人人爽秒播| 97超级碰碰碰精品色视频在线观看| 国产亚洲精品第一综合不卡| 亚洲成国产人片在线观看| 国产又爽黄色视频| 国产成人啪精品午夜网站|