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

    基于滑??刂频母呱璞蕊w行器協(xié)同制導(dǎo)方法

    2025-03-20 00:00:00郭博鐵鳴范文慧李傳旭
    關(guān)鍵詞:滑模控制參數(shù)優(yōu)化穩(wěn)定性

    摘 要:針對(duì)多個(gè)高升阻比飛行器在飛行中間段的時(shí)間協(xié)同問(wèn)題,提出一種基于滑模控制的高升阻比飛行器協(xié)同制導(dǎo)方法,建立多個(gè)高升阻比飛行器協(xié)同飛行力學(xué)模型,設(shè)計(jì)規(guī)劃飛行器在中間段制導(dǎo)過(guò)程中的飛行程序與針對(duì)飛行器側(cè)向機(jī)動(dòng)的滑??刂撇呗浴Mㄟ^(guò)設(shè)計(jì)標(biāo)稱軌跡,對(duì)控制參數(shù)進(jìn)行優(yōu)化,進(jìn)而生成各個(gè)飛行器不同初始條件所需的側(cè)向過(guò)載并通過(guò)側(cè)向過(guò)載得到所需的傾側(cè)角指令,以實(shí)現(xiàn)對(duì)飛行器飛行時(shí)間的控制,從而使多飛行器能夠同時(shí)到達(dá)目標(biāo)點(diǎn)上方設(shè)定范圍;考慮飛行器始末條件和狀態(tài)約束,使飛行器能夠滿足協(xié)同任務(wù)需求。利用李雅普諾夫穩(wěn)定性判據(jù)證明系統(tǒng)的穩(wěn)定性以及滑模面非奇異性,仿真結(jié)果表明,該協(xié)同制導(dǎo)策略具備一定抗干擾性,能夠滿足異地非同步發(fā)射的多個(gè)飛行器協(xié)同制導(dǎo)需求。

    關(guān)鍵詞: 高升阻比飛行器; 協(xié)同制導(dǎo); 滑??刂? 參數(shù)優(yōu)化; 穩(wěn)定性

    中圖分類號(hào): V 448.232

    文獻(xiàn)標(biāo)志碼: ADOI:10.12305/j.issn.1001 506X.2025.02.24

    Cooperative guidance method of high lift to drag ratio aircraft

    based on sliding mode control

    GUO Bo1, TIE Ming1,*, FAN Wenhui2, LI Chuanxu1

    (1. Science and Technology on Space Physics Laboratory, Beijing 100076, China;

    2. Department of Automation, Tsinghua University, Beijing 100084, China)

    Abstract:Aiming at the problem of time coordination of multiple high lift to drag ratio aircrafts in mid flight, a cooperative guidance method of high lift to drag ratio aircraft based on sliding mode control is proposed. A cooperative flight mechanics model of multiple lift to drag ratio aircrafts is established, and the flight program of the aircraft in mid flight guidance and the sliding mode control strategy for lateral maneuvering of the aircraft are designed and planned. By designing nominal trajectories, the control parameters are optimized to generate the lateral overload required by different initial conditions of each aircraft and obtain the required tilt angle command through the lateral overload, so as to control the flight time of the aircraft, so that multiple aircrafts can reach the set range above the target point at the same time, and consider the starting and ending conditions and state constraints of the aircraft, enabling the aircraft to meet collaborative mission requirements. Lyapunov stability criterion is used to prove the stability of the system and the non singularity of the sliding mode surface. Simulation results show that the cooperative guidance strategy has a certain anti interference property, and can meet the requirements of the cooperative guidance of multiple aircrafts launched asynchronously.

    Keywords:high lift to drag ratio aircraft; cooperative guidance; sliding mode control; parameter optimization; stability

    0 引 言

    高升阻比飛行器能夠在20~100 km的高層大氣飛行,相比一般的飛行器,該類飛行器飛行速度更快、飛行航程更遠(yuǎn)、飛行空域更廣,在一定時(shí)間內(nèi)改變飛行速度和飛行方向的能力更強(qiáng),因此研究該類飛行器具有重大意義1-2。然而,隨著工程需要,在實(shí)際場(chǎng)景中單個(gè)飛行器需逐漸勝任難度更高的任務(wù),多飛行器協(xié)同工作能力更強(qiáng),能夠勝任更高難度的任務(wù)。因此,多飛行器協(xié)同任務(wù)規(guī)劃在近些年來(lái)愈發(fā)受到關(guān)注3,需要考慮研究多飛行器協(xié)同的相關(guān)技術(shù)。

    對(duì)比單飛行器的任務(wù)規(guī)劃,多飛行器的協(xié)同規(guī)劃難度更高,需要考慮的約束更復(fù)雜,除了需要同時(shí)考慮多個(gè)飛行器的過(guò)載、動(dòng)壓、熱約束外4,還需要同時(shí)考慮在不同始末條件下多飛行器飛行時(shí)間的控制5

    實(shí)現(xiàn)多個(gè)飛行器針對(duì)同一目標(biāo)點(diǎn)的協(xié)同制導(dǎo)方式主要分為基于到達(dá)時(shí)間約束的飛行器獨(dú)立控制制導(dǎo)和基于信息交互的飛行器在線協(xié)同制導(dǎo)6-7。飛行器獨(dú)立時(shí)間控制制導(dǎo),即預(yù)先設(shè)定每個(gè)飛行器的飛行時(shí)間,在飛行過(guò)程中不進(jìn)行多個(gè)飛行器間的信息交互;基于信息交互的在線協(xié)同制導(dǎo)則利用多飛行器在飛行過(guò)程中的互相通信,實(shí)現(xiàn)時(shí)間協(xié)同8-9。而當(dāng)飛行器在大氣層高速飛行時(shí),飛行器周?chē)鷷?huì)形成等離子體鞘層,干擾通信10。因此,本文研究多飛行器的協(xié)同制導(dǎo),并采用更符合實(shí)際的獨(dú)立時(shí)間控制。

    針對(duì)多飛行器協(xié)同制導(dǎo)方面的研究,早期有些學(xué)者專門(mén)研究針對(duì)單個(gè)飛行器到達(dá)靜止目標(biāo)點(diǎn)的時(shí)間控制11-12。在此基礎(chǔ)上,又有更多協(xié)同制導(dǎo)方案被提出。文獻(xiàn)[13]首次將攻擊時(shí)間誤差反饋引入比例導(dǎo)引律,提出撞擊時(shí)間控制制導(dǎo)(impact time control guidance,ITCG),實(shí)現(xiàn)飛行器針對(duì)靜止目標(biāo)的協(xié)同。文獻(xiàn)[14]提出一種關(guān)于剩余時(shí)間和剩余距離的計(jì)算公式,通過(guò)預(yù)測(cè)剩余時(shí)間、剩余距離和末速度之間的關(guān)系建立能滿足到達(dá)時(shí)間約束的再入制導(dǎo)律。文獻(xiàn)[15]通過(guò)在橫向剖面使用預(yù)測(cè)校正制導(dǎo)調(diào)節(jié)時(shí)間以通過(guò)牛頓迭代法調(diào)節(jié)該剖面以滿足剩余距離和剩余時(shí)間的約束。王肖等人研究通過(guò)在高度-速度剖面內(nèi)設(shè)計(jì)參考軌跡,結(jié)合側(cè)向航向角走廊,通過(guò)在線數(shù)值預(yù)測(cè)剩余飛行航程和時(shí)間校正實(shí)現(xiàn)飛行器的時(shí)間協(xié)同制導(dǎo)16。周宏宇等人提出一種基于改進(jìn)粒子群優(yōu)化算法的多飛行器協(xié)同規(guī)劃算法,將強(qiáng)化學(xué)習(xí)方法同粒子群優(yōu)化算法結(jié)合,通過(guò)多次重構(gòu)飛行剖面最終完成協(xié)同制導(dǎo)17。You等人提出一種可以針對(duì)三維空間中多個(gè)飛行器的領(lǐng)航跟隨固定時(shí)間協(xié)同算法,使跟隨飛行器能夠通過(guò)算法,利用視線角度由領(lǐng)航飛行器帶領(lǐng)實(shí)現(xiàn)協(xié)同制導(dǎo)18。

    針對(duì)多飛行器協(xié)同控制,本文選取控制方式為滑??刂?。滑??刂苿?dòng)態(tài)響應(yīng)快,對(duì)參數(shù)攝動(dòng)不敏感,對(duì)外部的干擾具有較好的魯棒性,常被用來(lái)解決帶有碰撞角度約束的飛行器制導(dǎo)問(wèn)題19-20。文獻(xiàn)[21]提出一種考慮落角約束的滑模制導(dǎo)律,該方法通過(guò)改變分?jǐn)?shù)階參數(shù)實(shí)現(xiàn)對(duì)飛行器軌跡以及飛行時(shí)間的調(diào)節(jié)。文獻(xiàn)[22]和文獻(xiàn)[23]均為通過(guò)對(duì)法向過(guò)載使用滑??刂品椒ú⒃O(shè)定期望時(shí)間和期望角度,實(shí)現(xiàn)針對(duì)靜止目標(biāo)的協(xié)同制導(dǎo)。文獻(xiàn)[24]提出一種協(xié)同多個(gè)飛行器飛行時(shí)間和角度的三維自適應(yīng)滑模制導(dǎo)律,實(shí)現(xiàn)多飛行器在期望時(shí)刻的分布式協(xié)同制導(dǎo),并保證閉環(huán)系統(tǒng)的全局固定時(shí)間收斂特性。

    對(duì)于高升阻比飛行器來(lái)說(shuō),飛行器的動(dòng)力學(xué)模型都具有強(qiáng)非線性、強(qiáng)不確定性,其速度和加速度都是實(shí)時(shí)變化的。而且,相比末制導(dǎo)段,飛行器在中間段的飛行時(shí)間更長(zhǎng),距離更遠(yuǎn),有更足夠的空間來(lái)調(diào)節(jié)飛行時(shí)間25。因此,本文研究?jī)?nèi)容將進(jìn)一步考慮在氣動(dòng)力加速度隨大氣密度變化的情況下,多個(gè)飛行器在中間段的協(xié)同制導(dǎo),即多個(gè)飛行器從不同初始點(diǎn)開(kāi)始,依靠氣動(dòng)力通過(guò)滑模控制使飛行器在飛行時(shí)間相同的情況下到達(dá)同一目標(biāo)點(diǎn)上空,為末端制導(dǎo)做準(zhǔn)備。

    最后,通過(guò)仿真實(shí)驗(yàn)驗(yàn)證該方法的正確性和有效性。

    1 飛行器動(dòng)力學(xué)模型

    1.1 動(dòng)力學(xué)方程

    本文主要研究多飛行器在經(jīng)歷火箭助推段、自由段以及拉起段后在大氣層中無(wú)動(dòng)力飛行的中間段制導(dǎo)過(guò)程,假設(shè)第i枚飛行器的的動(dòng)力學(xué)模型可以表示如下:

    dvidt=-Dimi-gisin γi(1)

    dridt=visin γi(2)

    dNidt=vicos γisin iricos i(3)

    didt=vicos γicos iri(4)

    didt=Lisin σimivicos γi+Viricos γisinitan i(5)

    dγidt=Licos σimivi-givi-viricos γi(6)

    式中:下標(biāo)i代表第i枚飛行器;v為飛行器飛行速度;r為飛行器的地心距;N為飛行器所處經(jīng)度;為飛行器所處緯度;為飛行器航向角;γ為飛行器彈道傾角;D為飛行器受到的氣動(dòng)阻力;L為飛行器受到的氣動(dòng)升力;m為飛行器質(zhì)量;g為重力加速度大小;σ為飛行器傾側(cè)角。

    氣動(dòng)升力和氣動(dòng)阻力可表示如下:

    L=12ρv2SrefCL

    D=12ρv2SrefCD(7)

    式中:Sref為特征面積;CL和CD分別為升力系數(shù)和阻力系數(shù),可以分別由飛行器的攻角和馬赫數(shù)通過(guò)擬合計(jì)算得到26。

    h和ρ分別為飛行器所在高度和所在高度處的大氣密度:

    ρ=ρ0e-βh(8)

    h=r-R0(9)

    式中:ρ0為標(biāo)準(zhǔn)大氣密度;β為一常數(shù);R0為地球的平均半徑。

    取ρ0=1.225 kg/m3;β=1.406 4×10-4;R0=6 374.4 km。

    1.2 飛行約束條件

    每一枚飛行器在大氣層高速飛行的過(guò)程中,首先需要滿足氣動(dòng)熱約束、過(guò)載約束和動(dòng)壓約束:

    Qi=KQρiv3.15i≤Qmax(10)

    qi=0.5ρiv2i≤qmax(11)

    ni=L2i+D2imigi≤nmax(12)

    KQ=C1Rdρ0(g0R0)3.15(13)

    式中:Q為熱流密度;q為動(dòng)壓;n為總過(guò)載;KQ為熱流常數(shù);C1為形狀因子;Rd為前緣駐點(diǎn)處半徑;g0為地球表面加速度;Qmax為最大熱流密度;qmax為最大動(dòng)壓;nmax為最大過(guò)載。

    除此之外,還要求飛行器在滿足到達(dá)終端位置約束、終端高度約束、總航程約束的同時(shí),還需要保證多個(gè)飛行器到達(dá)的時(shí)間相同。

    設(shè)Srest,i為第i枚飛行器的剩余待飛航程,有以下關(guān)系式27

    Srest,i=R0·arccos[cos icos ′cos(θi-θ′)+

    sin isin ′](14)

    式中:i,θi分別表示第i枚飛行器當(dāng)前位置的緯度和經(jīng)度;′,θ′分別表示終點(diǎn)位置的緯度和經(jīng)度。

    設(shè)tf,i表示第i枚飛行器的到達(dá)時(shí)間,tdf為各個(gè)飛行器的期望到達(dá)時(shí)間,需要滿足關(guān)系:

    tf,1=tf,2=…=tf,i=tdf(15)

    2 多飛行器協(xié)同制導(dǎo)設(shè)計(jì)

    2.1 滑模面設(shè)計(jì)

    針對(duì)靜止目標(biāo)的制導(dǎo)幾何可以繪制在如圖1所示的參考系中。λ為飛行器相對(duì)目標(biāo)的視線角;γ為飛行器的航向角;θ表示航向角和視線角的差值,即飛行器前置角,前置角可以描述飛行器飛行方位角與目標(biāo)點(diǎn)的偏差,由此可以寫(xiě)成如下關(guān)系式:

    θ=λ-γ(16)

    vr=vcos θ(17)

    vλ=vsin θ(18)

    還可以由此得到視線角和航向角和速度分量vλ以及側(cè)向過(guò)載aM之間的關(guān)系:

    λ·=vλr(19)

    γ·=aMv(20)

    對(duì)于飛行器無(wú)動(dòng)力飛行的中段制導(dǎo)過(guò)程,飛行器的目標(biāo)是到達(dá)預(yù)定位置并與末制導(dǎo)段進(jìn)行交班。交班區(qū)域位置設(shè)計(jì)為在以目標(biāo)點(diǎn)經(jīng)緯度為圓心的圓形范圍內(nèi),當(dāng)飛行器到達(dá)該圓形范圍內(nèi)即完成任務(wù),飛行器還需在飛行末時(shí)刻滿足高度約束。

    針對(duì)多飛行器協(xié)同制導(dǎo),還需滿足多飛行器到達(dá)交班點(diǎn)的時(shí)間相同。因此,本文采用的是預(yù)測(cè)-校正制導(dǎo)方法,通過(guò)剩余飛行距離預(yù)估時(shí)間對(duì)飛行控制指令進(jìn)行實(shí)時(shí)校正以實(shí)現(xiàn)對(duì)飛行時(shí)間的控制28。

    假設(shè)交班點(diǎn)與終點(diǎn)的所剩距離為S;飛行器當(dāng)前位置距終點(diǎn)位置為S;飛行器到交班點(diǎn)所需的預(yù)估剩余時(shí)間為tgo;預(yù)估到達(dá)時(shí)刻為tf;當(dāng)前時(shí)刻為t;則可以得到如下關(guān)系式29

    tf=t+tgo(21)

    Srest=S-S(22)

    tgo=Srestv1+0.52N-1θ2(23)

    式中:N為比例導(dǎo)引系數(shù),一般設(shè)N=3??梢杂纱嗽O(shè)計(jì)滑模面為

    s=tf-tdf=t+tgo-tdf=tgo-tdgo(24)

    式中:tdf表示期望到達(dá)時(shí)間;tdgo表示期望剩余時(shí)間?;C鎠最終趨于穩(wěn)定。

    2.2 滑模制導(dǎo)律設(shè)計(jì)

    令r=Srest,通過(guò)對(duì)滑模面求導(dǎo),設(shè)計(jì)滑模制導(dǎo)律函數(shù)30

    s·=tgo·-tdgo·=r·v1+0.52N-1θ2-rθv(2N-1)λ·+

    1+rθv2(2N-1)am(25)

    基于式(25),可以利用李雅普諾夫非線性控制理論設(shè)計(jì)滑模控制的趨近律。

    aM=aeqM+aMconM+aswM(26)

    式中:aM為控制飛行器所需要的側(cè)向過(guò)載;aeqM為利用s·=0設(shè)計(jì)的等效控制分量,aeqM可以表示為

    aeqM=v2(2N-1)r·cos θ-1θ+0.5v2rθcos θ+VMλ·(27)

    同時(shí),為了使滑模面最終趨于零并保證系統(tǒng)滿足李雅普諾夫穩(wěn)定性條件,設(shè)計(jì)aswM和aMconM,可以寫(xiě)成如下關(guān)系式:

    aswM=-M(psgn(θ)+1)sgn(s)(28)

    aMconM=-h(huán)(θ)θks(29)

    在這里,可以引入一連續(xù)正函數(shù):

    h(θ)=sgn(θ)θ,|θ|lt;ε1

    1-ε1ε2-ε1|θ|+ε1(ε2-1)ε2-ε1, ε1≤|θ|≤ε2

    1, 其他(30)

    式中:k,s,M,p,ε1,ε2均為正常數(shù)且需要滿足pgt;1。

    接下來(lái),對(duì)系統(tǒng)的穩(wěn)定性進(jìn)行證明,證明過(guò)程可參考文獻(xiàn)[30],考慮李雅普諾夫函數(shù)為

    V=12s2(31)

    上述滑模面和趨近律的設(shè)計(jì),可以保證李雅普諾夫函數(shù)是正定的,接下來(lái)只需要證明其導(dǎo)數(shù)是負(fù)定,即可保證系統(tǒng)是穩(wěn)定的。

    對(duì)李雅普諾夫函數(shù)求導(dǎo)可得

    V·=ss·=-rh(θ)v2(2N-1)ks2-Mrv2(2N-1)(p|θ|+θ)|s|(32)

    分別取式(32)中兩項(xiàng)單獨(dú)討論,令

    A=rh(θ)v2(2N-1)(33)

    B=Mrv2(2N-1)(p|θ|+θ)|s|(34)

    由于h(θ)的函數(shù)值域均為正,k,M,p為正常數(shù)且pgt;1,所以當(dāng)θ≠0時(shí),A和B都是正定的。這意味著當(dāng)θ=0時(shí)不能滿足趨近律s=0,因此需要證明θ=0不是一個(gè)穩(wěn)定吸引子。

    下面考慮θ=0的情況,由式(16)求導(dǎo)得到:

    θ·=λ·-γ·(35)

    由式(18)~式(20)和式(26)可知,當(dāng)θ=0時(shí),能夠滿足λ·=0,aeqM=0,aMconM=0,由此可以得到:

    γ·=aswMv=-Mvsgn(s)(36)

    將式(36)和式(35)聯(lián)立得到:

    θ·=aswMv=Mvsgn(s)(37)

    因此,當(dāng)θ=0且s≠0時(shí),θ·≠0,系統(tǒng)并沒(méi)有趨于穩(wěn)定,所以θ=0并非一個(gè)穩(wěn)定性因子。僅當(dāng)s=0時(shí),θ=0才是一個(gè)吸引子。因此,由式(32)可以證明系統(tǒng)滿足李雅普諾夫穩(wěn)定性條件,接下來(lái)討論制導(dǎo)指令產(chǎn)生奇異性的問(wèn)題。理論上,當(dāng)r為零或無(wú)限接近于零時(shí)制導(dǎo)指令會(huì)產(chǎn)生奇異性現(xiàn)象,然而在實(shí)際飛行過(guò)程中,一般會(huì)將飛行結(jié)束條件設(shè)置為進(jìn)入S*范圍內(nèi)。此時(shí),Srest負(fù)定,Srest為零的狀態(tài)只是到達(dá)結(jié)束條件下的某一個(gè)短暫時(shí)刻,即不會(huì)在滑模面趨于穩(wěn)態(tài)時(shí)連續(xù)為零。因此,在式(26)中的滑模制導(dǎo)律不會(huì)導(dǎo)致奇異產(chǎn)生。同時(shí),另一個(gè)奇異條件是θ等于零或無(wú)限接近于零,而在式(27)中,令:

    T=cos θ-1θ(38)

    依據(jù)洛必達(dá)法則,當(dāng)θ無(wú)限趨于零的時(shí)候,T也是趨于零的,而不會(huì)變成無(wú)窮大項(xiàng),所以θ不會(huì)導(dǎo)致奇異情況產(chǎn)生。

    同樣,可以證明所需過(guò)載aM的另外兩個(gè)分量aMconM和aswM不會(huì)出現(xiàn)奇異,本文設(shè)計(jì)的滑模制導(dǎo)律不會(huì)出現(xiàn)奇異問(wèn)題。

    另外,為了抑制控制系統(tǒng)在運(yùn)行過(guò)程中產(chǎn)生的抖振現(xiàn)象,現(xiàn)將不連續(xù)函數(shù)sgn(x)改寫(xiě)為一連續(xù)函數(shù)sgmf(x):

    sgmf(x)=211+e-αz-12, αgt;0(39)

    2.3 參數(shù)優(yōu)化求解

    由于飛行器飛行速度、加速度都是實(shí)時(shí)變化的,所以為了防止飛行過(guò)程中傾側(cè)角變化出現(xiàn)抖振、能量不足等原因?qū)е聼o(wú)法在限定時(shí)間到達(dá)終端位置,可以對(duì)參數(shù)M進(jìn)行尋優(yōu)計(jì)算?,F(xiàn)采用提前規(guī)劃標(biāo)稱軌跡的方式,利用遞推最小二乘法對(duì)合適的M值進(jìn)行計(jì)算,其最小二乘指標(biāo)為

    J=min(tdf-tf,i)(40)

    本文采用遞推最小二乘法迭代求解合適的M值,利用飛行時(shí)間反向計(jì)算所需要的M值,以下給出使用該方法的計(jì)算過(guò)程。

    設(shè)計(jì)變量:

    ξ=M(41)

    建立如下殘差方程:

    g(ξ)=tdf-tf(42)

    由最小二乘法的極值原理,可以得到如下的迭代算法:

    ξk+1=ξk-A-1kBk(43)

    式中:

    Ak=∑Ni=1gTiξkgiξk(44)

    Bk=∑Ni=1gTiξkgi(ξk)(45)

    當(dāng)估計(jì)值與上一步的差的范數(shù)滿足一定精度時(shí)就停止迭代。

    2.4 縱向與側(cè)向軌跡規(guī)劃

    得到側(cè)向過(guò)載指令aM后,再由aM求得所需的程序角指令,由于飛行器采用傾側(cè)轉(zhuǎn)彎(bank to turn,BTT)控制方法,側(cè)向控制指令均由傾側(cè)角σ產(chǎn)生,根據(jù)側(cè)向過(guò)載aM,升力L和飛行器質(zhì)量m得到傾側(cè)角指令為

    σ=arcsin(aMmL)(46)

    由傾側(cè)角指令產(chǎn)生的側(cè)向過(guò)載是有限的,所以側(cè)向過(guò)載aM不能無(wú)限制地施加。另外,飛行器需要維持穩(wěn)定,因此需要對(duì)傾側(cè)角限幅,故限定傾側(cè)角范圍為σ∈[-60°,60°]。

    除傾側(cè)角外,飛行器的程序角還包括攻角,為滿足約束條件,還需要對(duì)飛行器的攻角α進(jìn)行規(guī)劃:

    α=amax, v≥v1

    aL/D max+amax-aL/D maxv1-v2(v-v2), v2lt;vlt;v1

    aL/D max, v≤v2(47)

    式中:amax為飛行器攻角的最大值;αL/D max為飛行器最大升阻比情況下的攻角,可由升阻比函數(shù)CL/CD求出;v1為一給定的較大速度值;v2為一給定的較小速度值。

    因此,所研究的問(wèn)題可以進(jìn)一步轉(zhuǎn)化為通過(guò)利用航程-傾側(cè)角剖面以及速度-攻角剖面設(shè)計(jì)的飛行走廊。

    飛行器的軌跡通過(guò)四階五步龍格庫(kù)塔法積分得到,每0.1 s進(jìn)行一次積分,當(dāng)滿足終端位置約束時(shí)結(jié)束積分,假設(shè)終端位置為進(jìn)入距離目標(biāo)點(diǎn)半徑為S的圓區(qū)域內(nèi)時(shí)即為飛行器完成任務(wù),最優(yōu)指標(biāo)為每一枚飛行器的到達(dá)時(shí)間一致且均為期望到達(dá)時(shí)間tfd,即

    J=min(tfd-tf,i), i=1,2,…,N(48)

    式中:k為飛行器的數(shù)量。

    3 仿真分析

    3.1 仿真背景及參數(shù)設(shè)定

    假設(shè)共有3枚飛行器且該3枚飛行器均為同種類的飛行器,各個(gè)飛行器的質(zhì)量均為1 000 kg,特征面積均設(shè)計(jì)為0.51 m2,氣動(dòng)參數(shù)參考使用CAV H(common aero vehicle high)飛行器的氣動(dòng)參數(shù),該類飛行器的過(guò)程約束可以參考文獻(xiàn)[31]并設(shè)置熱流上限為Qmax=1 200 kw/m2,過(guò)載上限為Nmax=4 g,動(dòng)壓上限為q=400 kPa。

    假設(shè)3枚飛行器(i=1,2,3)在經(jīng)過(guò)前面階段的飛行后于同一時(shí)刻開(kāi)始中間段飛行,分別由不同初始位置出發(fā),其初始條件分別如表1所示,而飛行器的目標(biāo)點(diǎn)與終端約束情況如表2所示。

    由表1可知,設(shè)定的3枚飛行器初始條件不同,因此其剩余航程均不相同。3枚飛行器剩余航程分別為4 939.4 km,5 062.0 km,5 005.2 km。飛行器初始前置角均為零,因此在初速度、初始高度、初始航跡傾角都相近的情況下,如果不通過(guò)合適側(cè)向機(jī)動(dòng)改變飛行時(shí)間,3枚飛行器最終到達(dá)終點(diǎn)的時(shí)間必不相同。

    3.2 多飛行器協(xié)同制導(dǎo)仿真分析

    對(duì)飛行器制導(dǎo)律的參數(shù)指標(biāo)進(jìn)行設(shè)計(jì),設(shè)k=0.001,p=1.3,ε1=0.001,ε2=0.015。

    對(duì)于期望到達(dá)時(shí)間tdf的設(shè)定有一定局限范圍,既不能過(guò)大也不能過(guò)小。因?yàn)榧纫WC其大于剩余航程最長(zhǎng)的飛行器在無(wú)側(cè)向機(jī)動(dòng)情況下的飛行時(shí)間,還要保證不能因飛行時(shí)間過(guò)長(zhǎng)導(dǎo)致能量提前消耗殆盡以致其無(wú)法到達(dá)終點(diǎn)。經(jīng)計(jì)算,能夠滿足3枚飛行器均可到達(dá)的期望時(shí)間可取范圍為tdf∈[1 366,1 404]s,因此設(shè)置tdf=1 386 s。通過(guò)第2節(jié)的最小二乘法計(jì)算出3枚飛行器的控制參數(shù)M值分別如表3所示。

    利用所設(shè)計(jì)的指標(biāo)進(jìn)行仿真,得到3枚飛行器在不同初始條件下的協(xié)同制導(dǎo)仿真結(jié)果。圖2展示了規(guī)劃出3枚飛行器的航程-高度飛行剖面,3枚飛行器在經(jīng)過(guò)不斷的跳躍飛行最終能夠完成所需的全部航程。圖3展示了3枚飛行器經(jīng)緯度變化的軌跡,可以看出飛行器由不同初始位置出發(fā),在經(jīng)過(guò)不斷側(cè)向機(jī)動(dòng)后最終能夠到達(dá)終點(diǎn)區(qū)域。

    圖4和圖5展示的是3枚飛行器飛行過(guò)程中高度和速度的變化。飛行器高度在不斷的來(lái)回波動(dòng)中緩慢下降至終端約束的高度范圍內(nèi),而飛行器速度隨著能量的消耗不斷減小。

    圖6展示了3枚飛行器飛行過(guò)程中動(dòng)壓、過(guò)載、熱流的變化情況,這些變化量在飛行過(guò)程中始終符合過(guò)程約束,不會(huì)超過(guò)限定值。

    圖7展示了3枚飛行器由制導(dǎo)指令產(chǎn)生的傾側(cè)角變化過(guò)程,為了改變飛行時(shí)間而調(diào)整飛行方向,傾側(cè)角需要根據(jù)制導(dǎo)指令在限定范圍內(nèi)實(shí)時(shí)調(diào)整。

    圖8展示了3枚飛行器在飛行過(guò)程中攻角的變化情況,為維持飛行高度以及實(shí)現(xiàn)更大航程,飛行器在飛行的過(guò)程中需要調(diào)整攻角,使攻角隨速度的衰減不斷增大。

    圖9和圖10分別展示了3枚飛行器的航向角變化和前置角變化,航向角和前置角在不斷地側(cè)向調(diào)整中,逐漸恢復(fù)到初始時(shí)的角度。

    3枚飛行器初始的航向角就是初始點(diǎn)與目標(biāo)點(diǎn)之間的方位角,因此在飛行器不進(jìn)行任何側(cè)向飛行的情況下,飛行器剛好能以最短時(shí)間準(zhǔn)確飛達(dá)目標(biāo)點(diǎn)上空。同理,初始前置角均接近于零,這意味著飛行器初始飛行方向是準(zhǔn)確朝著目標(biāo)點(diǎn)的,但是由于飛行器需要依靠側(cè)向機(jī)動(dòng)調(diào)整飛行時(shí)間,所以要進(jìn)行飛行方向的調(diào)整,而側(cè)向機(jī)動(dòng)的指令則由計(jì)算得到的傾側(cè)角實(shí)現(xiàn)。在完成對(duì)時(shí)間的調(diào)整后,前置角最終會(huì)逐漸恢復(fù)到零。航向角同樣如此,在完成對(duì)飛行時(shí)間的調(diào)整后也會(huì)逐漸恢復(fù)到初始航向角大小,進(jìn)而繼續(xù)準(zhǔn)確朝著目標(biāo)點(diǎn)進(jìn)發(fā)。

    最終,3枚飛行器均能在設(shè)定時(shí)間到達(dá)指定目標(biāo)區(qū)域內(nèi),其終端狀態(tài)如表4所示??梢?jiàn),飛行器均到達(dá)了距離目標(biāo)點(diǎn)140 km范圍處,3枚飛行器與期望到達(dá)時(shí)間的偏差分別為0.5 s、0.5 s、1 s??梢哉J(rèn)為在誤差范圍內(nèi),飛行器協(xié)同制導(dǎo)的任務(wù)是成功的。

    3.3 對(duì)比分析

    本文研究方法的思路主要是通過(guò)側(cè)向機(jī)動(dòng)以達(dá)到改變飛行器飛行時(shí)間的目的。而文獻(xiàn)[21]提出的方法主要針對(duì)基于縱向落角約束的預(yù)測(cè)-校正制導(dǎo)研究。因此,對(duì)比該方法的仿真如圖11~圖16所示。

    從結(jié)果看,文獻(xiàn)[21]方法通過(guò)預(yù)估的末時(shí)刻落角反饋制導(dǎo)指令,改變縱向飛行剖面。使用該方法同樣能夠達(dá)到協(xié)同制導(dǎo)的目的,并使飛行器在約束范圍內(nèi)到達(dá)距離目標(biāo)點(diǎn)140 km范圍內(nèi),完成交班任務(wù),最終的到達(dá)時(shí)間與期望到達(dá)時(shí)間的偏差分別為0.7 s、0.1 s、0.5 s,可以認(rèn)為在誤差范圍內(nèi)。

    但是,該方法通過(guò)改變制導(dǎo)律中的分?jǐn)?shù)階參數(shù),僅僅只改變飛行器的縱向軌跡,沒(méi)有通過(guò)橫側(cè)向制導(dǎo)改變飛行器的橫側(cè)向軌跡,因此飛行器只可始終沿初始航向角飛行。而飛行器在實(shí)際飛行中難免受到干擾,不能在無(wú)側(cè)向機(jī)動(dòng)控制的情況下始終保持初始方位角,最終可能造成飛行方向的偏差。

    相比文獻(xiàn)[21]中的方法,本文研究方法更優(yōu)異的方面在于能夠通過(guò)直接對(duì)飛行器的橫側(cè)向軌跡進(jìn)行控制約束,使飛行器在飛行過(guò)程中糾正到達(dá)時(shí)間偏差的同時(shí),還可以對(duì)飛行方向進(jìn)行修正,最終令飛行器始終能夠在制導(dǎo)指令下向目標(biāo)位置飛行。

    除此之外,利用與本文研究方法同樣的思路,使用傳統(tǒng)的比例-積分-微分(proportion integration differentiation,PID)控制器法設(shè)置一組對(duì)照實(shí)驗(yàn)。該方法與本文研究思路類似,同樣利用預(yù)估剩余時(shí)間tgo和期望剩余時(shí)間tdgo設(shè)計(jì)關(guān)于側(cè)向過(guò)載的PID制導(dǎo)律:

    aM=KP(tgo-tdgo)+KI∫tf0(tgo-tdgo)dt+

    KD(tgo·-tdgo·)(49)

    式中:KP、KI、KD分別為PID控制參數(shù)。利用與前面所述處理滑??刂茀?shù)相同的方法,固定KI、KD,通過(guò)最小二乘迭代計(jì)算調(diào)節(jié)參數(shù)KP。

    圖17~圖19分別展示了在傳統(tǒng)PID控制方法下的3枚飛行器航程-高度飛行剖面、過(guò)程約束以及傾側(cè)角指令的變化。由圖17可以看出,3枚飛行器最終在飛行終點(diǎn)處均能夠到達(dá)指定范圍內(nèi)的高度。由圖18可以看出,3枚飛行器均能在過(guò)程約束范圍內(nèi)飛行。但是,通過(guò)圖19得到,3枚飛行器的傾側(cè)角指令在飛行最后時(shí)段發(fā)生了在上下幅值之間的高頻率來(lái)回抖振,這意味著需要控制傾側(cè)角變化的飛行器舵發(fā)生高頻率來(lái)回?cái)[動(dòng),顯然是不符合實(shí)際的。另外,在使用該傳統(tǒng)PID控制方法得到的結(jié)果中,3枚飛行器到達(dá)時(shí)刻分別為1 388.2 s、1 387.4 s、1 386.0 s,與期望到達(dá)時(shí)間的最大偏差達(dá)到了2.2 s,結(jié)果指標(biāo)略遜于本文所述方法。因此,通過(guò)仿真實(shí)驗(yàn)可知,本文使用方法相比于傳統(tǒng)PID控制方法要更優(yōu)。

    3.4 蒙特卡羅仿真分析

    在實(shí)際飛行過(guò)程中,飛行器會(huì)受到各種不確定性因素的干擾,所以需要對(duì)飛行器在受到外界擾動(dòng)情況下的飛行條件進(jìn)行仿真驗(yàn)證。因此,為了進(jìn)一步驗(yàn)證該方法的抗干擾性,選取一組數(shù)據(jù)利用蒙特卡羅打靶法進(jìn)行仿真。該組數(shù)據(jù)為滿足正態(tài)分布條件的初始偏差和外界擾動(dòng),如表5所示,偏差限為3σ。

    接下來(lái)以飛行器3為例,驗(yàn)證本文研究方法中飛行器在微小擾動(dòng)情況下的實(shí)際效果。通過(guò)100次蒙特卡羅打靶仿真,得到如下結(jié)果。

    圖20為飛行器在擾動(dòng)條件下高度隨時(shí)間變化情況,飛行器飛行過(guò)程中高度受到影響,產(chǎn)生微小變化,但是在終點(diǎn)仍能符合預(yù)定高度約束范圍。圖21為飛行器到達(dá)時(shí)間與目標(biāo)到達(dá)時(shí)間tdf偏差的密度分布,在滿足正態(tài)分布的偏差條件下,飛行器到達(dá)時(shí)間偏差同樣符合正態(tài)分布,與預(yù)期情況一致。但是正態(tài)分布的均值不為零,原因是最初由最小二乘法得到的最優(yōu)控制參數(shù)M計(jì)算出的飛行器到達(dá)時(shí)間本就有微小偏差,而此正態(tài)分布的均值即為最初到達(dá)時(shí)間與目標(biāo)到達(dá)時(shí)間tdf的差。但是,該正態(tài)分布的時(shí)間偏差的最大值為1.7 s,仍然在可以接受的誤差范圍內(nèi)。因此,可以認(rèn)為該多飛行器協(xié)同制導(dǎo)方法能克服干擾因素的影響。

    圖22展示了飛行器到達(dá)時(shí)間偏差與其對(duì)應(yīng)的末速度,直觀展示了飛行器末速度隨時(shí)間偏差的變化,這些偏差均在可接受范圍內(nèi)。

    4 結(jié) 論

    針對(duì)多高升阻比飛行器的協(xié)同制導(dǎo)技術(shù)是未來(lái)復(fù)雜環(huán)境與任務(wù)亟需發(fā)展的重點(diǎn)研究方向。本文提出一種通過(guò)預(yù)估剩余時(shí)間和預(yù)估剩余航程的多飛行器協(xié)同制導(dǎo)方法,該方法可以在保證制導(dǎo)準(zhǔn)確的同時(shí)調(diào)節(jié)飛行時(shí)間,實(shí)現(xiàn)多飛行器協(xié)同制導(dǎo)。

    通過(guò)設(shè)計(jì)基于滑模制導(dǎo)律并利用最小二乘迭代法計(jì)算控制參數(shù),規(guī)劃了飛行器飛行所需攻角和傾側(cè)角指令。通過(guò)攻角指令維持飛行高度,使飛行器能在到達(dá)終端狀態(tài)時(shí)滿足高度約束,通過(guò)傾側(cè)角指令使飛行器能有效調(diào)節(jié)飛行時(shí)間,解決復(fù)雜條件下對(duì)飛行器飛行時(shí)間控制的問(wèn)題,使飛行器能在多約束條件下實(shí)現(xiàn)協(xié)同制導(dǎo)。最后,通過(guò)仿真實(shí)驗(yàn)驗(yàn)證了該方法的可行性。仿真結(jié)果表明,多個(gè)飛行器的到達(dá)時(shí)間與預(yù)期到達(dá)時(shí)間的差距均能控制在1 s內(nèi)。對(duì)比其他方法和思路,本文研究方法具有更優(yōu)異的表現(xiàn)。并且,通過(guò)蒙特卡羅仿真實(shí)驗(yàn)可得,即使在有初始條件偏差和過(guò)程擾動(dòng)的情況下,飛行器的時(shí)間偏差仍然不超過(guò)2 s,說(shuō)明本方法具有一定抗干擾性,為高升阻比飛行器協(xié)同制導(dǎo)問(wèn)題的解決提供了一種有效的技術(shù)手段。

    參考文獻(xiàn)

    [1]DING Y B, YUE X K, CHEN G S, et al. Review of control and guidance technology on hypersonic vehicle[J]. Chinese Journal of Aeronautics, 2022, 35(7): 1-18.

    [2]WU T C, WANG H L, LIU Y H, et al. Learning based interfered fluid avoidance guidance for hypersonic reentry vehicles with multiple constraints[J]. ISA Transactions, 2023, 139(8): 291-307.

    [3]JIA X, WU S T, WEN Y M, et al. A distributed decision method for missiles autonomous formation based on potential game[J]. Journal of Systems Engineering and Electronics, 2019, 30(4): 738-748.

    [4]LIANG Z X, REN Z, LI Q D. Evolved atmospheric entry corridor with safety factor[J]. Acta Astronautica, 2018, 143(2): 82-91.

    [5]ZHANG D, LIU L, WANG Y J. On line reentry guidance algorithm with both path and no fly zone constraints[J]. Acta Astronautica, 2015, 117(12): 243-253.

    [6]YOU H, CHANG X L, ZHAO J F, et al. Three dimensional impact angle constrained cooperative guidance strategy against maneuvering target[J]. ISA Transactions, 2023, 138(7): 262-280.

    [7]ZHAO J, ZHOU R, DONG Z N. Three dimensional cooperative guidance laws against stationary and maneuvering targets[J]. Chinese Journal of Aeronautics, 2015,28(4): 1104-1120.

    [8]SONG J H, SONG S M, XU S L. Three dimensional cooperative guidance law for multiple missiles with finite time convergence[J]. Aerospace Science and Technology, 2017, 67(8): 193-205.

    [9]WANG X H, LU X. Three dimensional impact angle constrai ned distributed guidance law design for cooperative attacks[J]. ISA Transactions, 2018, 73(2): 79-90.

    [10]LI J F, ZHOU J X, YAO J F, et al. Experimental observations of communication in blackout,topological waveguiding and Dirac zero index property in plasma sheath[J]. Nanophotonics, 2023, 12(10): 1847-1856.

    [11]LEE J I, JEON I S, TAHK M J. Guidance law to control impact time and angle[J]. IEEE Trans.on Aerospace and Electronic Systems, 2007, 43(1): 301-310.

    [12]KIM T H, LEE C H, JEON I S, et al. Augmented polynomial guidance with impact time and angle constraints[J]. IEEE Trans.on Aerospace and Electronic Systems, 2013, 49(4): 2806-2817.

    [13]JEON I S, LEE J I, TAHK M J. Impact time control guidance law for anti ship missiles[J]. IEEE Trans.on Control Systems Technology, 2006, 14(2): 260-266.

    [14]GUO Y H, LI X, ZHANG H J, et al. Entry guidance with terminal time control based on quasi equilibrium glide condition[J]. IEEE Trans.on Aerospace and Electronic Systems, 2019, 56(2): 887-896.

    [15]LI Z H, HE B, WANG M H, et al. Time coordination entry guidance for multi hypersonic vehicles[J]. Aerospace Science and Technology, 2019, 89(6): 123-135.

    [16]王肖, 郭杰, 唐勝景, 等. 基于解析剖面的時(shí)間協(xié)同再入制導(dǎo)[J]. 航空學(xué)報(bào), 2019, 40(3): 239-250.

    WANG X, GUO J, TANG S J, et al. Time cooperative entry guidance based on analytical profile[J]. Acta Aeronautica et Astronautica Sinica, 2019, 40(3): 239-250.

    [17]周宏宇, 王小剛, 單永志, 等. 基于改進(jìn)粒子群算法的飛行器協(xié)同軌跡規(guī)劃[J]. 自動(dòng)化學(xué)報(bào), 2022, 48(11): 2670-2676.

    ZHOU H Y, WANG X G, SHAN Y Z, et al. Synergistic path planning for multiple vehicles based on an improved particle swarm optimization method[J]. Acta Automatica Sinica, 2022, 48(11): 2670-2676.

    [18]YOU H, CHANG X L, ZHAO J F, et al. Three dimensional impact angle constrained fixed time cooperative guidance algorithm with adjustable impact time[J]. Aerospace Science and Technology, 2023, 141: 108574.

    [19]ZHANG W J, FU S N, LI W, et al. An impact angle constraint integral sliding mode guidance law for maneuvering targets interception[J]. Journal of Systems Engineering and Electronics, 2020, 31(1): 168-184.

    [20]ZHAO F J, YOU H. New three dimensional secondorder sliding mode guidance law with impact angle constraints[J]. The Aeronautical Journal, 2020, 124(1273): 368-384.

    [21]盛永智, 甘佳豪, 張成新. 彈道可調(diào)的落角約束分?jǐn)?shù)階滑模制導(dǎo)律設(shè)計(jì)[J]. 航空學(xué)報(bào), 2023, 44(7): 182-195.

    SHENG Y Z, GAN J H, ZHANG C X. Fractional order sliding mode guidance law design with trajectory adjustable and terminal angular constraint[J]. Acta Aeronautica et Astronautica Sinica, 2023, 44(7): 182-195.

    [22]LI W, WEN Q Q, HE L, et al. Three dimensional impact angle constrained distributed cooperative guidance law for anti ship missiles[J]. Journal of Systems Engineering and Electro nics, 2021, 32(2): 447-459.

    [23]LYU T, LI C J, GUO Y N, et al.Three dimensional finite time cooperative guidance for multiple missiles without radial velocity measurements[J]. Chinese Journal of Aeronautics, 2019, 32(5): 1294-1304.

    [24]王雨辰, 王偉, 林時(shí)堯, 等. 考慮攻擊時(shí)間及空間角度約束的三維自適應(yīng)滑模協(xié)同制導(dǎo)律設(shè)計(jì)[J]. 兵工學(xué)報(bào), 2023, 44(9): 2778-2790.

    WANG Y C, WANG W, LIN S Y, et al. Three dimensional adaptive sliding mode cooperative guidance law with impact time and angle constraints[J]. Acta Armamentarii, 2023, 44(9): 2778-2790.

    [25]BAO C Y, WANG P, TANG G J. Integrated method of gui dance,control and morphing for hypersonic morphing vehicle in glide phase[J]. Chinese Journal of Aeronautics, 2021, 34(5): 535-553.

    [26]SHENG Y Z, ZHANG Z, XIA L. Fractional order sliding mode control based guidance law with impact angle constraint[J]. Nonlinear Dynamics, 2021, 106(1): 425-444.

    [27]姜鵬, 郭棟, 韓亮, 等. 多飛行器再入段時(shí)間協(xié)同彈道規(guī)劃方法[J]. 航空學(xué)報(bào), 2020, 41(S1): 171-183.

    JIANG P, GUO D, HAN L, et al. Trajectory optimization for cooperative reentry of multiple hypersonic glide vehicle[J]. Acta Aeronautica et Astronautica Sinica, 2020, 41(S1): 171-183.

    [28]LIANG Z X, LV C, ZHU S. Lateral entry guidance with terminal time constraint[J]. IEEE Trans.on Aerospace and Electronic System, 2023, 59(3): 2544-2553.

    [29]JEON I S, LEE J I, TAHK M J. Homing guidance law for cooperative attack of multiple missiles[J]. Journal of Guidance, Control, and Dynamics, 2010, 33(1): 275-280.

    [30]CHO D, KIM H J, TAHK M J. Nonsingular sliding mode guidance for impact time control[J]. Journal of Guidance, Control, and Dynamics, 2016, 39(1): 61-68.

    [31]徐慧, 蔡光斌, 崔亞龍, 等. 高超聲速滑翔飛行器再入軌跡優(yōu)化[J]. 哈爾濱工業(yè)大學(xué)學(xué)報(bào), 2023, 55(4): 44-55.

    XU H, CAI G B, CUI Y L, et al. Reentry trajectory optimization method of hypersonic glide vehicle[J]. Journal of Harbin Institute of Technology, 2023,55(4): 44-55.

    作者簡(jiǎn)介

    郭 博(1999—),男,碩士研究生,主要研究方向?yàn)轱w行器制導(dǎo)與控制、飛行器智能決策規(guī)劃與仿真。

    鐵 鳴(1976—),男,研究員,博士,主要研究方向?yàn)橹悄茱w行器設(shè)計(jì)、復(fù)雜系統(tǒng)建模與仿真。

    范文慧(1966—),男,教授,博士,主要研究方向?yàn)閺?fù)雜系統(tǒng)建模與仿真。

    李傳旭(1997—),男,工程師,碩士,主要研究方向?yàn)轱w行器制導(dǎo)與控制、系統(tǒng)仿真。

    猜你喜歡
    滑??刂?/a>參數(shù)優(yōu)化穩(wěn)定性
    非線性中立型變延遲微分方程的長(zhǎng)時(shí)間穩(wěn)定性
    基于干擾觀測(cè)器的PID滑模變結(jié)構(gòu)控制
    基于神經(jīng)網(wǎng)絡(luò)的動(dòng)力電池組焊接參數(shù)優(yōu)化研究
    基于多算法的ROV遠(yuǎn)程協(xié)同控制器設(shè)計(jì)與實(shí)現(xiàn)
    半動(dòng)力系統(tǒng)中閉集的穩(wěn)定性和極限集映射的連續(xù)性
    研究LTE與WCDMA系統(tǒng)間小區(qū)互操作與參數(shù)優(yōu)化
    基于磁流變技術(shù)的汽車(chē)發(fā)動(dòng)機(jī)隔振系統(tǒng)的參數(shù)優(yōu)化
    科技視界(2016年23期)2016-11-04 08:17:36
    改進(jìn)的動(dòng)態(tài)面船舶航跡跟蹤控制研究
    上向進(jìn)路式尾砂膠結(jié)充填采礦法采場(chǎng)結(jié)構(gòu)參數(shù)優(yōu)化研究
    中國(guó)科技博覽(2016年5期)2016-04-23 05:45:45
    国产成人精品婷婷| 国产欧美亚洲国产| 老汉色av国产亚洲站长工具| 亚洲 欧美一区二区三区| 少妇的丰满在线观看| 国产精品人妻久久久影院| 亚洲色图 男人天堂 中文字幕| 国产有黄有色有爽视频| 一级毛片黄色毛片免费观看视频| 五月开心婷婷网| 人人妻人人爽人人添夜夜欢视频| 亚洲国产欧美网| 亚洲精品国产av蜜桃| 老司机影院成人| 久久综合国产亚洲精品| 在线天堂中文资源库| 老汉色∧v一级毛片| 欧美精品国产亚洲| 两个人免费观看高清视频| 人成视频在线观看免费观看| 18禁裸乳无遮挡动漫免费视频| 最近手机中文字幕大全| av视频免费观看在线观看| 十八禁网站网址无遮挡| 1024视频免费在线观看| 天堂俺去俺来也www色官网| 国产在视频线精品| 免费观看av网站的网址| 免费高清在线观看日韩| 日韩,欧美,国产一区二区三区| 久久久久久久久久人人人人人人| 国产福利在线免费观看视频| av在线观看视频网站免费| 免费播放大片免费观看视频在线观看| 中文字幕色久视频| 国产视频首页在线观看| 啦啦啦视频在线资源免费观看| 中文字幕精品免费在线观看视频| 国产精品免费视频内射| 黄色 视频免费看| 在线观看人妻少妇| 久久亚洲国产成人精品v| 十八禁网站网址无遮挡| 免费女性裸体啪啪无遮挡网站| 日韩av在线免费看完整版不卡| 看非洲黑人一级黄片| 国产精品亚洲av一区麻豆 | 日韩精品免费视频一区二区三区| 一区福利在线观看| 欧美变态另类bdsm刘玥| 十八禁高潮呻吟视频| 久久久久网色| 飞空精品影院首页| 黄网站色视频无遮挡免费观看| 国语对白做爰xxxⅹ性视频网站| 精品少妇一区二区三区视频日本电影 | 精品少妇内射三级| 国产熟女午夜一区二区三区| 日韩大片免费观看网站| 日韩精品有码人妻一区| 两性夫妻黄色片| 一本色道久久久久久精品综合| 日韩,欧美,国产一区二区三区| 人人妻人人爽人人添夜夜欢视频| 国产日韩欧美亚洲二区| 午夜福利视频在线观看免费| 少妇熟女欧美另类| 男女午夜视频在线观看| 日韩一区二区视频免费看| 亚洲欧美一区二区三区久久| 国产片特级美女逼逼视频| 五月天丁香电影| 久久热在线av| 午夜福利影视在线免费观看| 免费观看性生交大片5| 男男h啪啪无遮挡| 午夜福利视频精品| 国产精品秋霞免费鲁丝片| 国产一区二区在线观看av| 久久久久久久久久久免费av| 日韩一本色道免费dvd| a 毛片基地| 久久精品国产综合久久久| 男的添女的下面高潮视频| 国产精品一区二区在线观看99| 国产免费又黄又爽又色| 成人18禁高潮啪啪吃奶动态图| 亚洲国产日韩一区二区| 丝袜在线中文字幕| 高清黄色对白视频在线免费看| 99国产精品免费福利视频| 成年女人毛片免费观看观看9 | 爱豆传媒免费全集在线观看| 久久久久久伊人网av| 国产精品 国内视频| 欧美bdsm另类| 最近最新中文字幕大全免费视频 | 午夜老司机福利剧场| 蜜桃国产av成人99| 成年av动漫网址| 美女视频免费永久观看网站| 涩涩av久久男人的天堂| 免费黄频网站在线观看国产| 成人午夜精彩视频在线观看| 日韩,欧美,国产一区二区三区| 日韩中文字幕欧美一区二区 | 国产精品成人在线| 成年女人毛片免费观看观看9 | 亚洲精品国产av成人精品| 亚洲四区av| 精品少妇一区二区三区视频日本电影 | 高清视频免费观看一区二区| 国产白丝娇喘喷水9色精品| 精品人妻在线不人妻| 色网站视频免费| 亚洲av欧美aⅴ国产| 免费大片黄手机在线观看| 这个男人来自地球电影免费观看 | 欧美日韩一区二区视频在线观看视频在线| 成人毛片60女人毛片免费| 黄片小视频在线播放| 色网站视频免费| 十八禁高潮呻吟视频| 久久av网站| 捣出白浆h1v1| 伊人亚洲综合成人网| 亚洲精品日韩在线中文字幕| 边亲边吃奶的免费视频| 永久网站在线| 大香蕉久久网| 免费人妻精品一区二区三区视频| 人人妻人人澡人人爽人人夜夜| 日韩中字成人| 熟女少妇亚洲综合色aaa.| 国产成人91sexporn| 国产日韩一区二区三区精品不卡| 久久国内精品自在自线图片| 曰老女人黄片| 狂野欧美激情性bbbbbb| 黑人巨大精品欧美一区二区蜜桃| 久热这里只有精品99| 99热全是精品| 一个人免费看片子| 精品人妻在线不人妻| 老汉色av国产亚洲站长工具| 欧美最新免费一区二区三区| 日本av免费视频播放| 亚洲欧洲国产日韩| 大码成人一级视频| 大片免费播放器 马上看| 纯流量卡能插随身wifi吗| 日本欧美国产在线视频| 日本91视频免费播放| 97精品久久久久久久久久精品| 精品99又大又爽又粗少妇毛片| 国产亚洲精品第一综合不卡| 国产精品蜜桃在线观看| 国产不卡av网站在线观看| 亚洲精品国产av蜜桃| 老汉色av国产亚洲站长工具| 久久av网站| 在线天堂最新版资源| 在线天堂最新版资源| 精品一区二区免费观看| 国产精品一区二区在线不卡| 亚洲成国产人片在线观看| 高清不卡的av网站| 日日爽夜夜爽网站| 最近2019中文字幕mv第一页| 丁香六月天网| 一个人免费看片子| 少妇被粗大的猛进出69影院| 人人妻人人爽人人添夜夜欢视频| 成年美女黄网站色视频大全免费| 黄色一级大片看看| 欧美亚洲 丝袜 人妻 在线| 国产女主播在线喷水免费视频网站| 国产精品成人在线| 国产男女内射视频| 久久99精品国语久久久| 欧美少妇被猛烈插入视频| 色视频在线一区二区三区| av女优亚洲男人天堂| 久久久国产欧美日韩av| 性色avwww在线观看| 成年av动漫网址| 国产白丝娇喘喷水9色精品| 男女免费视频国产| 精品久久蜜臀av无| 国产成人免费观看mmmm| 丰满少妇做爰视频| 国产在线一区二区三区精| 美女午夜性视频免费| 狠狠婷婷综合久久久久久88av| 国产精品偷伦视频观看了| 欧美xxⅹ黑人| 国产精品秋霞免费鲁丝片| 午夜福利在线观看免费完整高清在| 日韩精品有码人妻一区| 亚洲av电影在线观看一区二区三区| 丰满迷人的少妇在线观看| 嫩草影院入口| 99热全是精品| 久久久久网色| 在线观看一区二区三区激情| av在线播放精品| 婷婷成人精品国产| 亚洲,欧美精品.| 午夜日本视频在线| 亚洲av国产av综合av卡| 十分钟在线观看高清视频www| 国精品久久久久久国模美| 亚洲欧洲国产日韩| av电影中文网址| 日韩制服丝袜自拍偷拍| 久久久久久久久久久免费av| 成人国产麻豆网| 亚洲色图综合在线观看| 秋霞在线观看毛片| 国产探花极品一区二区| 国产成人精品无人区| 日产精品乱码卡一卡2卡三| 国产人伦9x9x在线观看 | 黄频高清免费视频| 亚洲国产欧美在线一区| 黄色毛片三级朝国网站| 日韩中文字幕欧美一区二区 | 制服诱惑二区| 成人国语在线视频| 亚洲精品国产av成人精品| 亚洲精品av麻豆狂野| 熟妇人妻不卡中文字幕| 欧美日韩国产mv在线观看视频| 伊人久久大香线蕉亚洲五| 免费黄网站久久成人精品| 久久久久国产网址| 午夜福利一区二区在线看| 汤姆久久久久久久影院中文字幕| 91国产中文字幕| 黄色怎么调成土黄色| 亚洲欧美中文字幕日韩二区| 亚洲欧美一区二区三区久久| 午夜精品国产一区二区电影| 亚洲,欧美精品.| 两性夫妻黄色片| 91aial.com中文字幕在线观看| 国精品久久久久久国模美| 日韩欧美精品免费久久| 欧美av亚洲av综合av国产av | 99热全是精品| 精品人妻熟女毛片av久久网站| 亚洲精品自拍成人| 夫妻午夜视频| 免费在线观看视频国产中文字幕亚洲 | 999精品在线视频| 亚洲男人天堂网一区| 一级片免费观看大全| 午夜激情久久久久久久| 侵犯人妻中文字幕一二三四区| 一二三四在线观看免费中文在| 尾随美女入室| 汤姆久久久久久久影院中文字幕| 亚洲美女黄色视频免费看| 日日爽夜夜爽网站| 国产女主播在线喷水免费视频网站| 亚洲精品,欧美精品| 国产成人精品婷婷| 欧美国产精品va在线观看不卡| 成年动漫av网址| 国产熟女欧美一区二区| 国产男人的电影天堂91| 午夜日韩欧美国产| 国产在线一区二区三区精| 日本黄色日本黄色录像| 99国产精品免费福利视频| 日韩精品有码人妻一区| 丰满少妇做爰视频| 国语对白做爰xxxⅹ性视频网站| 999久久久国产精品视频| 欧美+日韩+精品| 一区二区三区精品91| 久久99热这里只频精品6学生| 国产精品女同一区二区软件| 国产精品久久久久久久久免| 国产xxxxx性猛交| 午夜日本视频在线| 啦啦啦在线观看免费高清www| 一级黄片播放器| 亚洲精品美女久久av网站| 中文天堂在线官网| 99热国产这里只有精品6| 女性被躁到高潮视频| 久久久久久人妻| 晚上一个人看的免费电影| 精品国产露脸久久av麻豆| 中文天堂在线官网| 午夜免费鲁丝| 99久久人妻综合| 亚洲欧美精品自产自拍| 少妇的丰满在线观看| 亚洲五月色婷婷综合| 欧美老熟妇乱子伦牲交| 国产精品国产三级专区第一集| 亚洲欧美成人综合另类久久久| 中文欧美无线码| 日日摸夜夜添夜夜爱| 久久久久久人妻| 伦精品一区二区三区| 久久久精品国产亚洲av高清涩受| 精品第一国产精品| 国产深夜福利视频在线观看| 国产一区二区激情短视频 | 99精国产麻豆久久婷婷| 七月丁香在线播放| 亚洲欧美精品综合一区二区三区 | 黄片小视频在线播放| 日韩一本色道免费dvd| 一二三四中文在线观看免费高清| 日韩制服丝袜自拍偷拍| 国产一区二区激情短视频 | 伦理电影免费视频| 精品第一国产精品| 侵犯人妻中文字幕一二三四区| 在线观看人妻少妇| 老鸭窝网址在线观看| 国产日韩一区二区三区精品不卡| 自拍欧美九色日韩亚洲蝌蚪91| 国产无遮挡羞羞视频在线观看| 久久av网站| 极品人妻少妇av视频| 又黄又粗又硬又大视频| kizo精华| 婷婷色综合大香蕉| 亚洲精品国产一区二区精华液| 中文乱码字字幕精品一区二区三区| 三上悠亚av全集在线观看| 色网站视频免费| 精品国产一区二区久久| 久久精品夜色国产| 亚洲欧美一区二区三区久久| 欧美在线黄色| 欧美日韩视频精品一区| 亚洲精品久久午夜乱码| av在线老鸭窝| 老女人水多毛片| 成人手机av| 久久久久久久精品精品| xxxhd国产人妻xxx| 欧美在线黄色| 免费播放大片免费观看视频在线观看| 有码 亚洲区| 久久精品国产亚洲av天美| 人人妻人人澡人人爽人人夜夜| 精品国产乱码久久久久久男人| 成人漫画全彩无遮挡| 日韩一区二区三区影片| 成年人免费黄色播放视频| 国产在线视频一区二区| av网站在线播放免费| 女人被躁到高潮嗷嗷叫费观| 久久精品国产鲁丝片午夜精品| 熟女电影av网| 啦啦啦在线观看免费高清www| 在线亚洲精品国产二区图片欧美| 丝瓜视频免费看黄片| 亚洲精品视频女| 久久97久久精品| 午夜福利乱码中文字幕| 深夜精品福利| 丰满乱子伦码专区| 久久人人爽人人片av| 午夜免费观看性视频| 亚洲美女视频黄频| 最近2019中文字幕mv第一页| 韩国精品一区二区三区| 超色免费av| 午夜福利视频精品| 亚洲图色成人| 中文天堂在线官网| 美女福利国产在线| 老女人水多毛片| 1024视频免费在线观看| 久久国产亚洲av麻豆专区| 人妻 亚洲 视频| 国产白丝娇喘喷水9色精品| √禁漫天堂资源中文www| 大片免费播放器 马上看| 99国产精品免费福利视频| 波多野结衣一区麻豆| 青春草亚洲视频在线观看| 一级爰片在线观看| 一级黄片播放器| 亚洲三级黄色毛片| 中文字幕人妻熟女乱码| 亚洲欧美中文字幕日韩二区| 欧美av亚洲av综合av国产av | 青草久久国产| 国产 精品1| 中文字幕精品免费在线观看视频| 国产精品国产av在线观看| 日韩大片免费观看网站| 黄片播放在线免费| 国产色婷婷99| 亚洲欧美成人综合另类久久久| 色婷婷久久久亚洲欧美| 久久精品人人爽人人爽视色| xxx大片免费视频| 美女福利国产在线| 精品亚洲成国产av| av福利片在线| 亚洲国产精品一区三区| 一区福利在线观看| 成人手机av| 久久精品夜色国产| 宅男免费午夜| 精品久久蜜臀av无| 自线自在国产av| freevideosex欧美| 国产精品女同一区二区软件| 免费少妇av软件| 日韩制服骚丝袜av| 三级国产精品片| 国产精品二区激情视频| 一级毛片 在线播放| 中文字幕制服av| 亚洲国产看品久久| 少妇人妻 视频| 免费黄色在线免费观看| 香蕉精品网在线| 大香蕉久久成人网| 性色avwww在线观看| 国产精品.久久久| 90打野战视频偷拍视频| 欧美日韩视频高清一区二区三区二| 婷婷色av中文字幕| 成人毛片a级毛片在线播放| 美女福利国产在线| 亚洲精品美女久久av网站| 国产亚洲最大av| 性色avwww在线观看| 五月开心婷婷网| 国产精品一区二区在线不卡| 久久久久国产一级毛片高清牌| 久久精品国产亚洲av天美| 一区在线观看完整版| 成年人午夜在线观看视频| 肉色欧美久久久久久久蜜桃| 色视频在线一区二区三区| 国产精品久久久av美女十八| 欧美日韩成人在线一区二区| 在线免费观看不下载黄p国产| 深夜精品福利| av天堂久久9| 男人添女人高潮全过程视频| 中文字幕色久视频| 夜夜骑夜夜射夜夜干| 亚洲国产日韩一区二区| 国产精品女同一区二区软件| 国产精品国产av在线观看| 五月开心婷婷网| 久久97久久精品| 国产在线视频一区二区| 亚洲av福利一区| 中文乱码字字幕精品一区二区三区| 国产熟女欧美一区二区| 国语对白做爰xxxⅹ性视频网站| 国产爽快片一区二区三区| 两性夫妻黄色片| 男女边吃奶边做爰视频| 久久久久久久久免费视频了| 波多野结衣av一区二区av| 一区在线观看完整版| 老汉色∧v一级毛片| 亚洲av日韩在线播放| 日本免费在线观看一区| 性色avwww在线观看| 老司机影院毛片| 国产免费一区二区三区四区乱码| 国产在线一区二区三区精| 免费播放大片免费观看视频在线观看| 中文字幕制服av| 91午夜精品亚洲一区二区三区| 亚洲综合精品二区| 亚洲成人av在线免费| 亚洲婷婷狠狠爱综合网| 人妻 亚洲 视频| 中文精品一卡2卡3卡4更新| 国产色婷婷99| 国产男人的电影天堂91| 丰满饥渴人妻一区二区三| 欧美亚洲 丝袜 人妻 在线| 美国免费a级毛片| 香蕉丝袜av| 亚洲av中文av极速乱| 一边摸一边做爽爽视频免费| 男女边摸边吃奶| 午夜激情av网站| a级毛片在线看网站| 美女xxoo啪啪120秒动态图| 熟女av电影| av天堂久久9| 国产成人午夜福利电影在线观看| 久久这里只有精品19| 国产乱人偷精品视频| 亚洲av男天堂| 大片免费播放器 马上看| 精品酒店卫生间| 午夜福利视频精品| 久久久a久久爽久久v久久| 国产精品人妻久久久影院| 美女午夜性视频免费| 又粗又硬又长又爽又黄的视频| 美女国产高潮福利片在线看| 人成视频在线观看免费观看| 一区二区三区四区激情视频| 在线观看免费视频网站a站| 亚洲精品国产色婷婷电影| 久热久热在线精品观看| 亚洲天堂av无毛| a级毛片黄视频| 国产精品 国内视频| 国产亚洲午夜精品一区二区久久| 80岁老熟妇乱子伦牲交| 久久久久人妻精品一区果冻| 免费黄色在线免费观看| 丝袜人妻中文字幕| 日韩一区二区三区影片| 国产成人精品久久久久久| 老司机影院成人| 少妇 在线观看| 2022亚洲国产成人精品| 99久国产av精品国产电影| 久久久亚洲精品成人影院| 亚洲av电影在线观看一区二区三区| 老汉色∧v一级毛片| av卡一久久| 成人亚洲欧美一区二区av| 在线观看国产h片| 老汉色av国产亚洲站长工具| 国产免费视频播放在线视频| 亚洲精品美女久久久久99蜜臀 | 在线观看国产h片| 亚洲国产日韩一区二区| av国产精品久久久久影院| 女性生殖器流出的白浆| 性色avwww在线观看| 寂寞人妻少妇视频99o| 国产一区二区 视频在线| 不卡视频在线观看欧美| av卡一久久| 黄色一级大片看看| 午夜91福利影院| 久久精品国产综合久久久| 大话2 男鬼变身卡| 性色avwww在线观看| 最近最新中文字幕大全免费视频 | 青春草视频在线免费观看| av片东京热男人的天堂| 美女中出高潮动态图| 激情五月婷婷亚洲| 伊人久久国产一区二区| 亚洲国产av新网站| 超碰97精品在线观看| 搡女人真爽免费视频火全软件| 纯流量卡能插随身wifi吗| 日韩伦理黄色片| 一级a爱视频在线免费观看| 亚洲av欧美aⅴ国产| 久久精品久久久久久久性| 日韩不卡一区二区三区视频在线| av国产久精品久网站免费入址| 亚洲一码二码三码区别大吗| 国产日韩欧美视频二区| 免费高清在线观看视频在线观看| 18禁动态无遮挡网站| 亚洲av中文av极速乱| 国产精品国产三级国产专区5o| 最新中文字幕久久久久| 波多野结衣av一区二区av| 国产精品欧美亚洲77777| 人人妻人人爽人人添夜夜欢视频| 丰满迷人的少妇在线观看| 在线精品无人区一区二区三| 热re99久久国产66热| 日本-黄色视频高清免费观看| 久久久久久久精品精品| 人体艺术视频欧美日本| 久久久久久久精品精品| 熟女电影av网| 嫩草影院入口| 女性生殖器流出的白浆| 老司机影院成人| 91国产中文字幕| 国产精品免费视频内射| 美女主播在线视频| 不卡av一区二区三区| 亚洲精品aⅴ在线观看| 黄网站色视频无遮挡免费观看| 视频区图区小说| 少妇熟女欧美另类| 黄片无遮挡物在线观看| 精品亚洲乱码少妇综合久久| 少妇人妻 视频| 亚洲精品久久成人aⅴ小说| 丰满乱子伦码专区| 亚洲伊人色综图| 黑人猛操日本美女一级片| 欧美激情极品国产一区二区三区| 夜夜骑夜夜射夜夜干| 欧美亚洲 丝袜 人妻 在线| 王馨瑶露胸无遮挡在线观看| 永久网站在线| 国产精品国产三级国产专区5o| 精品人妻偷拍中文字幕| 亚洲精品成人av观看孕妇|