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

    基于有限元法的水下航行器地磁異常模擬研究

    2024-11-21 00:00:00趙高陽劉勇朱平杰向冰周洪娟
    關(guān)鍵詞:磁場模型

    摘" 要:

    為改進(jìn)水下航行器磁異常探測技術(shù)中磁異常實(shí)測數(shù)據(jù)難獲得、現(xiàn)有磁體模型計(jì)算精度不高等問題,基于有限元仿真探討建立無需區(qū)分近、遠(yuǎn)場的高精度混合磁體模型方法。應(yīng)用有限元數(shù)值方法仿真復(fù)雜結(jié)構(gòu)水下航行器的空間磁異常分布,以磁場數(shù)值解為收斂目標(biāo),建立以偶極子陣元數(shù)目、位置及磁矩為參數(shù)的均勻磁化橢球體與偶極子陣列混合的航行器磁異常解析模型,采用非線性最小二乘算法求解模型系數(shù)。仿真結(jié)果表明,基于該模型得到磁異常計(jì)算值與全空間內(nèi)數(shù)值解的擬合度較高,測試平面平均誤差為3%。該模型在磁場延拓、高精度建模等方面可以進(jìn)一步應(yīng)用。關(guān)鍵詞:

    航空磁異常探測; 有限元仿真; 磁體模擬; 數(shù)值擬合

    中圖分類號(hào):

    O 411.3

    文獻(xiàn)標(biāo)志碼: A""" DOI:10.12305/j.issn.1001-506X.2024.07.03

    Simulation of geomagnetic anomaly of underwater vehicle based on

    finite element method

    ZHAO Gaoyang1, LIU Yong2, ZHU Pingjie2, XIANG Bing2, ZHOU Hongjuan1,*

    (1. School of Information Science and Engineering, Harbin Institute of Technology (Weihai), Weihai 264209, China;

    2. 722 Research Institute, China State Shipbuilding Corporation, Wuhan 430205, China)

    Abstract:

    To overcome the difficulties of obtaining the magnetic anomaly data of underwater vehicles and the low calculation precision of existing magnet simulation models, a high-precision hybrid analytical model suitable in both near-field and far-field based on finite element method is introduced. The spatial signature of magnetic anomaly of an underwater vehicle with complex structure is analyzed by the finite element method. Then, an analytical model combined with uniformly-magnetized spheroid and a dipole array is proposed with the numerical results as the convergence goal. This analytical model consists of several parameters, including the number of array elements, the positions of each element, and magnetic moments, which are solved by nonlinear least-square algorithm. The simulation results show that the magnetic anomaly calculated from this model has good consistencies with the numerical solution in full space, giving an average error of 3% at the measurement plane. This model could be potentially extended to applications like magnetic field continuation and high-precision modeling.

    Keywords:

    airborne magnetic anomaly detection; finite element method (FEM); magnet simulation; numerical fitting

    \=

    0" 引" 言

    水下航行器在日益重要的海洋開發(fā)中承擔(dān)著水下目標(biāo)、資源探測、氣象、海洋數(shù)據(jù)搜集等任務(wù)[1],這使得對水下航行器的探測具有深刻的民用和軍事意義。對水下大型航行器的探測來說,隨著非磁性或弱磁性材料的應(yīng)用與降噪技術(shù)的發(fā)展,傳統(tǒng)聲學(xué)探測受到前所未有的挑戰(zhàn),特別是在對聲信號(hào)具有巨大影響的淺海區(qū)域[23]。航空磁異常探測技術(shù)利用機(jī)載磁探儀對鐵磁性目標(biāo)進(jìn)行探測與定位,具有全天時(shí)和全地域的特點(diǎn)。作為一種非接觸式的被動(dòng)探測方法[45],航空磁異常探測具有識(shí)別能力強(qiáng)、定位精度高、隱蔽性良好等優(yōu)點(diǎn)[6],在一定程度上能彌補(bǔ)水聲探測的應(yīng)用限制。

    常用的水下航行器磁異常估算方法有數(shù)值法和磁體模擬解析法。數(shù)值法主要包括邊界元法、有限元法等[78],可以較好地解決不規(guī)則物體的磁場模擬問題,缺點(diǎn)是計(jì)算量大且實(shí)時(shí)性差,無法適應(yīng)實(shí)際應(yīng)用。磁體模擬解析法通過特定磁矩磁源產(chǎn)生的磁場來模擬目標(biāo)磁場[910]。常用的模型有單磁偶極子、磁偶極子陣列、均勻磁化的旋轉(zhuǎn)橢球體等[11]。目前,普遍采用的模擬方法是:在靠近航行器一倍艇長的距離內(nèi),航行器近似為單個(gè)旋轉(zhuǎn)橢球,其精度為75%~85%[12];當(dāng)距離遠(yuǎn)至可忽略航行器物理尺寸時(shí),航行器可簡化為單個(gè)磁偶極子或磁偶極子陣列模型,其中對磁偶極子陣列模型還涉及到多個(gè)磁偶極子的磁矩估計(jì)問題??紤]近、遠(yuǎn)場模型的切換問題,出現(xiàn)了磁偶極子和旋轉(zhuǎn)橢球體的組合磁體模型,常用的組合磁體擬合方法需要求解多維超定方程,存在穩(wěn)定性與方程系數(shù)矩陣復(fù)雜的問題,同時(shí)建模過程中偶極子數(shù)量及位置改變也會(huì)影響方程求解。

    在確定磁體模型磁矩值方面,目前多數(shù)研究基于實(shí)測數(shù)據(jù)或采用經(jīng)驗(yàn)值來確定,如一艘中型航行器(艇長為70~80 m)磁矩值為2×108 mA·m2[13]。采用單一經(jīng)驗(yàn)值的做法在理論研究上可以得到具有指導(dǎo)意義的定性結(jié)論,但在實(shí)際應(yīng)用中具有精度不高、靈活性不足的缺點(diǎn)。為獲得精確磁矩值,已有針對遺傳算法[14]、模擬退火算法[15]等智能算法的相關(guān)研究,但在討論中以固定的磁偶極子個(gè)數(shù)為前提并基于固定數(shù)目數(shù)據(jù)點(diǎn)討論建模誤差[16],不能達(dá)到全面優(yōu)化磁矩值的目的,模擬方法仍然缺乏靈活性。

    本文在已有的航行器磁異常磁體模擬的研究基礎(chǔ)上[17],基于COMSOL Multiphysics多物理場軟件(簡稱為COMSOL軟件)仿真航行器目標(biāo)磁異常分布,通過有限元仿真方法獲得數(shù)值解以解決難以獲取實(shí)測數(shù)據(jù)的困難,進(jìn)一步建立參數(shù)包括可變陣元數(shù)目的偶極子陣列、磁偶極子位置以及磁矩的旋轉(zhuǎn)橢球體與偶極子陣列混合磁體模型,以仿真數(shù)值解為收斂目標(biāo),實(shí)現(xiàn)具有良好穩(wěn)定性并且不區(qū)分近、遠(yuǎn)場探測的高精度磁體模擬方法的數(shù)值擬合研究。

    1" 混合磁體模型

    由磁疇磁化理論可知,組成航行器外殼體及內(nèi)部設(shè)備儀器的鐵磁性材料內(nèi)部自發(fā)磁化形成磁疇結(jié)構(gòu)[1819],在地磁場的作用下各磁疇的排列方向近似一致[20],形成感應(yīng)磁場,該感應(yīng)場由內(nèi)部磁矩方向相同的磁疇結(jié)構(gòu)的感應(yīng)場疊加而成[21],因此在磁體模擬法等效模型中,對各個(gè)磁體磁矩可選取相同的磁化方向。

    為提高航行器磁異常的估算精度,建立偶極子陣列、均勻磁化旋轉(zhuǎn)橢球體的混合磁體模型(以下簡稱混合模型),如圖1所示。

    在圖1中,以旋轉(zhuǎn)橢球體的中心為原點(diǎn),n元偶極子陣列等間距排布在旋轉(zhuǎn)橢球體長軸,橢球體長半軸為a,短半軸為b,與目標(biāo)的幾何尺寸相當(dāng)。

    混合磁體在探測點(diǎn)產(chǎn)生的磁感應(yīng)強(qiáng)度B表示為以下矩陣形式:

    B(r,n,M)=∑nj=0Bj=∑nj=0ajxajyajz

    bjxbjybjz

    cjxcjycjzMj(1)

    式中:r為空間探測點(diǎn)相對于原點(diǎn)或偶極子中心的位置矢量;n為偶極子陣列陣元個(gè)數(shù);Mj表示磁矩,其中j=0代表旋轉(zhuǎn)橢球體的磁矩矢量,1≤j≤n時(shí)表示第j個(gè)磁偶極子的磁矩矢量;Mjx、Mjy、Mjz分別表示沿x、y、z方向的磁矩分量。由旋轉(zhuǎn)橢球體和磁偶極子磁場的推導(dǎo),可得式(1)中相關(guān)系數(shù)[22]表示如下。

    (1) j=0

    ajx=-3μ04π12g2lnan+gan-g-ang2t,ajy=3μ0xy4πanb2nt,ajz=3μ0xz4πanb2nt

    bjx=ajy,bjy=-3μ08πanb2ng2-12g2lnan+gan-g-2any2b4nt,bjz=3μ0anxz4πanb2nt

    cjx=ajz,cjy=bjz,cjz=-3μ08πanb2ng2-12g2lnan+gan-g-2anz2b4nt

    式中:an=12(x2+y2+z2+g2+t);bn=an-g2t=(x2+y2+z2+g2)2-4g2x2;g=a2-b2。

    (2) 1≤j≤n

    ajx=14π3(x-xj)2r5j-1r3j,ajy=3(x-xj)y4πr5j,ajz=3(x-xj)z4πr5j

    bjx=ajy,bjy=14π3y2r5j-1r3j,bjz=3yz4πr5j

    cjx=ajz,cjy=bjz,cjz=14π3z2r5j-1r3j

    式中:rj=(x-xj)2+(y-yj)2+(z-zj)2;(xj,yj,zj)為第j個(gè)偶極子的位置坐標(biāo)。

    為提高擬合精度,由混合模型的磁感應(yīng)場分量Bx、By、Bz組成含待求系數(shù)的方程組:

    Bx=∑nj=0(ajxMjx+ajyMjy+ajzMjz)(2)

    By=∑nj=0(bjxMjx+bjyMjy+bjzMjz)(3)

    Bz=∑nj=0(cjxMjx+cjyMjy+cjzMjz)(4)

    由電磁場的唯一性定理可知,求解電磁場的定解問題時(shí),只要滿足給定的控制方程和邊界條件,其解必定唯一。因此,在求解中將探測點(diǎn)(x,y,z)坐標(biāo)值及對應(yīng)磁場值代入式(2)~式(4),即可建立含偶極子陣元數(shù)目、對應(yīng)磁矩參數(shù)、陣元位置的非線性方程組。

    2" 航行器磁異常的有限元仿真

    磁體模擬法需要部分磁場數(shù)據(jù),目前難以獲得實(shí)際海域中目標(biāo)磁異常實(shí)測數(shù)據(jù),且實(shí)測數(shù)據(jù)往往存在外界干擾、儀器噪聲甚至較多無效數(shù)據(jù),因此利用數(shù)值計(jì)算法獲得磁異常的仿真應(yīng)用具有較強(qiáng)的理論研究價(jià)值。文獻(xiàn)[23]中討論了以球體驗(yàn)證有限元軟件COMSOL仿真結(jié)果與理論值的一致性,為簡化運(yùn)算,將通過該軟件對目標(biāo)磁特性進(jìn)行仿真。有限元計(jì)算原理是求解如下偏微分方程。

    在無電流區(qū)域內(nèi),靜磁場的基本方程如下:

    Δ×H=0(5)

    定義標(biāo)量磁勢Vm的表達(dá)式為

    H=-ΔVm(6)

    設(shè)目標(biāo)殼體材料的相對磁導(dǎo)率為μr,利用磁通密度與磁場的本構(gòu)關(guān)系,可得到:

    Δ·B=0(7)

    B=μ0μrH(8)

    由此,得到標(biāo)量磁勢方程:

    -Δ·(μ0μrΔVm)=0(9)

    在本研究背景下,水下航行器外加磁場為地磁場Bt,待求解的微分方程為

    -Δ·(μ0μrΔVm+Bt)=0(10)

    按照艦長為100 m、艇直徑為10 m的中型航行器建立幾何模型,對殼體參數(shù)、材料相對磁導(dǎo)率、背景磁場等參數(shù)盡量按照真實(shí)參數(shù)設(shè)置。如圖2所示,構(gòu)建水下航行器模型并進(jìn)行網(wǎng)格劃分,完成計(jì)算。

    圖3以探測高度為500 m的平面目標(biāo)磁感應(yīng)場分布為例。受到目標(biāo)物理形狀的影響,無論沿橫向還是縱向,其感應(yīng)磁場的分布并不滿足嚴(yán)格對稱分布。在傳統(tǒng)磁體模擬方法中,通常采用單個(gè)均勻磁化橢球體或偶極子模型來進(jìn)行近似模擬,模擬的磁感應(yīng)場呈對稱分布,與實(shí)際磁異常分布不符,勢必導(dǎo)致在某些探測位置上磁場強(qiáng)度失真較大。傳統(tǒng)的單磁體模擬法對理論研究具有重要的意義,但不適用于對磁異常測量精度要求較高的實(shí)際應(yīng)用場合,因此構(gòu)建合適的混合磁體模型實(shí)現(xiàn)對目標(biāo)感應(yīng)磁場的高精度擬合非常必要。

    3" 混合磁體模型系數(shù)求解及驗(yàn)證

    3.1" 系數(shù)求解原理及過程

    本文采用非線性最小二乘算法列文伯格馬夸爾特(Levenberg-Marquarat, L-M)算法對待定系數(shù)進(jìn)行求解[24]。非線性最小二乘算法的原理是從某個(gè)初值開始按照搜索方向持續(xù)迭代計(jì)算,再通過新方向的調(diào)整,直至目標(biāo)函數(shù)收斂到極小值點(diǎn),較常用的算法有最速下降法[25]、牛頓法、高斯牛頓法等。上述算法在本文中均存在一定的使用缺陷:最速下降法在接近于最優(yōu)解時(shí)收斂速度較慢,因此需要較多迭代次數(shù);牛頓法需計(jì)算目標(biāo)函數(shù)的二階導(dǎo)矩陣,在數(shù)據(jù)規(guī)模較大時(shí)計(jì)算困難[26];高斯牛頓法因近似二階泰勒展開,只在展開點(diǎn)附近有較好結(jié)果,在迭代中可能出現(xiàn)奇異或病態(tài)矩陣,或者因局部近似無法實(shí)現(xiàn)收斂[2728]。作為改進(jìn),L-M算法為增量添加信任域,該區(qū)域范圍依近似模型與實(shí)際函數(shù)的差異自適應(yīng)調(diào)整[29]。

    目標(biāo)函數(shù)定義為磁體模型與實(shí)際函數(shù)間誤差平方和R:

    R=∑Ni=1[f2(ri)-f1(ri,n,M)]2(11)

    式中:i=0,1,2,…,N表示用以求解的探測點(diǎn);ri為探測距離;n為偶極子陣元數(shù)目;M為對應(yīng)磁矩系數(shù);f1表示含待求參數(shù)的磁體模型;f2表示實(shí)際值函數(shù),在本文中為COMSOL軟件磁場仿真值。該算法通過多次迭代計(jì)算,控制目標(biāo)函數(shù)最小以獲得待求參數(shù),主要包括以下步驟。

    步驟 1" f1待求參數(shù)p,初始點(diǎn)p0,終止控制常數(shù)ε0=f2-f1(pk),參數(shù)k=0,λ0=10-3,v可以是任何大于1的數(shù)。

    步驟 2" Jacobi矩陣為Jk,構(gòu)造增量正規(guī)方程N(yùn)-kδk=JTkεk,其中N-k=JTkJk+λkI。

    步驟 3" 求解增量正規(guī)方程得到δk。

    步驟 3.1" 若f2-f1(pk+δk)lt;εk,令pk+1=pk+εk,當(dāng)δklt;ε時(shí)輸出結(jié)果;否則迭代繼續(xù),并令λk+1=λk/v。

    步驟 3.2" 如果f2-f1(pk+δk)≥εk,令λk+1=λk·v,重新解正規(guī)方程得到δk,返回步驟1。

    如圖4所示,水下航行器外殼形狀可能不完全規(guī)則,參照其形狀特征建立適當(dāng)尺寸橢球體并取其長軸方向上的中心點(diǎn)為原點(diǎn)??臻g探測點(diǎn)Qi對應(yīng)磁感應(yīng)場Bi=(Bix, Biy, Biz)。

    圖4中,rij=(rijx, rijy, rijz,)為探測點(diǎn)Qi與磁偶極子陣列第j個(gè)陣元相對距離。

    建立混合模型含參方程組如下,利用最小二乘算法求解待求磁矩參數(shù),并建立橢球體與磁偶極子陣列的混合磁體模型。

    Bix=∑nj=0(ajxMjx+ajyMjy+ajzMjz)

    Biy=∑nj=0(bjxMjx+bjyMjy+bjzMjz)

    Biz=∑nj=0(cjxMjx+cjyMjy+cjzMjz)(12)

    3.2" 擬合實(shí)驗(yàn)結(jié)果與分析

    下面將通過實(shí)驗(yàn)驗(yàn)證擬合模型。通過L-M算法求解混合模型中的待定參數(shù),在擬合過程中可根據(jù)應(yīng)用條件需求

    設(shè)定陣元數(shù)目范圍、對最佳陣元數(shù)進(jìn)行取整求解,選擇與探測點(diǎn)磁場值誤差最小的L-M算法結(jié)果,算法初值參考經(jīng)驗(yàn)值。后文進(jìn)行了多組實(shí)驗(yàn)討論取點(diǎn)方式、陣元數(shù)目等參數(shù)影響,增大模型驗(yàn)證范圍。假定水下航行器艇長為100 m,以建立橢球體與五元偶極子陣列混合模型為例,基于數(shù)值計(jì)算結(jié)果隨機(jī)選取探測點(diǎn),分別以完全規(guī)律性、一定規(guī)律性、完全隨機(jī)性3種方式選取探測點(diǎn),討論選取探測點(diǎn)的規(guī)律性對最終模型擬合精度的影響。

    (1) 擬合使用符合完全規(guī)律性探測點(diǎn)

    在固定探測高度平面上選取規(guī)律分布探測點(diǎn)如圖5所示,在500 m探測高度平面上,以均勻間隔在與航行器運(yùn)行方向夾角為45°的航線上選取探測點(diǎn)Qi=(x,y,500),x,y∈(-500,500),i=0,1,2,…,100。將對應(yīng)位置坐標(biāo)與磁場值代入式(12),并利用L-M算法對磁矩參數(shù)進(jìn)行擬合。

    在500 m探測高度平面上,對探測平面磁場分布的數(shù)值仿真值與混合模型計(jì)算值等效效果進(jìn)行評(píng)估,分別繪制磁感應(yīng)場與3個(gè)分量,如圖6所示。計(jì)算磁場計(jì)算值與仿真值間峰值差:x分量差值為0.050 nT,y分量差值為0.039 nT,z分量差值為0.066 nT,磁場總量差值為0.046 nT,通過繪制場強(qiáng)等值線分布圖可知磁場分布擬合結(jié)果較接近。

    (2) 擬合使用符合一定規(guī)律性探測點(diǎn)

    在固定探測高度平面上隨機(jī)選取相同數(shù)目探測點(diǎn)如圖7所示,探測點(diǎn)Qi=(x, y,500),x,y∈(-500,500),i=0,1,2,…,100;將對應(yīng)位置坐標(biāo)與磁場值代入式(11)進(jìn)行擬合。

    進(jìn)行后續(xù)的磁場分布計(jì)算,同樣選取500 m探測高度測平面,對該取點(diǎn)方式計(jì)算得到混合模型磁場分布計(jì)算值,與數(shù)值仿真值等效效果進(jìn)行評(píng)估,分別繪制磁感應(yīng)場與3個(gè)分量,如圖8所示。計(jì)算磁場計(jì)算值與仿真值間峰值差:x分量差值為0.089 nT,y分量差值為0.029 nT,z分量差值為0.052 nT,磁場總量差值為0.004 nT,雖然場值分布有差異但差值較小。

    (3) 擬合使用符合完全隨機(jī)性探測點(diǎn)

    更一般地,在選取探測點(diǎn)時(shí)不限制探測高度,在探測空間中完全隨機(jī)選擇相同數(shù)目探測點(diǎn),Qi=(x,y,z),x,y∈(-500,500),z∈(400,1 200),i=0,1,2,…,100。在本文的討論中,探測空間由磁探儀作用距離進(jìn)行定義:國外領(lǐng)先技術(shù)將有效作用范圍提升到1 200 m左右;有磁探儀據(jù)稱有效作用距離大約為400~600 m,因此將有效探測距離設(shè)定為航行器上方400~1 200 m距離。在該空間內(nèi)隨機(jī)選取探測點(diǎn),如圖9所示,利用探測點(diǎn)坐標(biāo)信息與磁場數(shù)值進(jìn)行擬合。

    在模型計(jì)算值與數(shù)值仿真值的對比中,同樣選取500 m探測高度測平面,分別繪制磁感應(yīng)場與3個(gè)分量等效效果,如圖10所示。計(jì)算磁場計(jì)算值與仿真值間峰值差:x分量差值為0.038 nT,y分量差值為0.033 nT,z分量差值為0.057 nT,磁場總量差值為0.016 nT。

    在混合模型建模時(shí)分別用3種選取探測點(diǎn)方法計(jì)算待求參數(shù),3種不同模型對目標(biāo)空間磁場分布計(jì)算值均較接近仿真值。圖11以平均絕對百分比誤差(mean absolute percentage error, MAPE)量化目標(biāo)磁場分布計(jì)算值與仿真值間誤差,在計(jì)算誤差的過程中擴(kuò)大探測平面數(shù)據(jù)密度,驗(yàn)證點(diǎn)坐標(biāo)范圍為x,y∈(-1 000, 1 000),z∈(400,1 200)。

    分析不同取點(diǎn)方式,建立模型計(jì)算誤差,3種模型計(jì)算誤差均較穩(wěn)定,應(yīng)用在不同模型使用條件下可適當(dāng)調(diào)整擬合取點(diǎn)方式:采用完全規(guī)律與一定規(guī)律取點(diǎn)方法,擬合所用探測點(diǎn)來自同一高度平面,在該平面及附近高度平面內(nèi)誤差較小。采用完全隨機(jī)取點(diǎn)方法擬合所用探測點(diǎn)數(shù)據(jù)范圍更大,增強(qiáng)了模型的魯棒性,計(jì)算誤差比較穩(wěn)定,特別是應(yīng)用在較遠(yuǎn)探測距離時(shí)使模型誤差降低,對后續(xù)數(shù)值模型的應(yīng)用與研究有一定的啟發(fā)作用。

    3.3" 擬合結(jié)果應(yīng)用

    聯(lián)合有限元仿真與磁體模擬方法獲得高精度的數(shù)值擬合模型,理論上可以通過航行器目標(biāo)引起地磁異常分布,判斷大致位置。

    以圖6為例進(jìn)行分析,航行器形狀因?yàn)榫哂袑ΨQ、規(guī)則結(jié)構(gòu),磁感應(yīng)場及3個(gè)分量具有明顯對稱性;x分量具有很強(qiáng)的區(qū)域特征,在航行器位置處出現(xiàn)1個(gè)負(fù)峰值;y分量有4個(gè)較突出磁場峰值,其中2個(gè)正峰值分別出現(xiàn)在目標(biāo)的左前側(cè)和右后側(cè)、2個(gè)負(fù)峰值出現(xiàn)在左后側(cè)和右前側(cè),且右側(cè)峰值較左邊更大,右側(cè)峰值位置以航行器位置為中心沿y軸對稱;磁異常z方向分量峰值最大,對目標(biāo)地磁異常總場影響最大,在探測平面上各有1個(gè)磁場正峰值與負(fù)峰值,正峰值位置出現(xiàn)在航行器運(yùn)行方向的前側(cè),與磁異??偭康?/p>

    峰值位置大致相同,負(fù)峰出現(xiàn)在后側(cè),目標(biāo)大概位置在磁場負(fù)峰值中心與正峰值中心中間距離處。按照以上的分布特征,由磁場分布及磁場分量特征可確定航行器位于x∈(-80, 50)、y=0,與實(shí)際仿真條件大致相同,驗(yàn)證了該方法在確定目標(biāo)位置上的可行性。

    4" 混合模型建模精度分析

    4.1" 驗(yàn)證模型精度改善

    本文在建模方法上對單個(gè)磁體模擬方法進(jìn)行改進(jìn),經(jīng)前文對比與討論,實(shí)驗(yàn)中建立的五元偶極子陣列混合模型實(shí)現(xiàn)較好擬合。為進(jìn)一步對建模精度進(jìn)行分析,在相同的探測條件下對以下磁體模型等效性進(jìn)行評(píng)價(jià)。

    Ⅰ模型:均勻磁化的單個(gè)旋轉(zhuǎn)橢球體模型。

    Ⅱ模型:橢球體與五元偶極子陣列混合模型。

    Ⅲ模型:橢球體與九元偶極子陣列混合模型。

    求解磁體模型中的待求系數(shù)。Ⅰ模型中未知磁矩值參考常用經(jīng)驗(yàn)公式,Ⅱ、Ⅲ模型通過L-M算法求解參數(shù),在建模時(shí)保證模型結(jié)構(gòu)尺寸等相同。建模完成后,在500 m探測高度測平面上分別計(jì)算3個(gè)模型磁感應(yīng)場及分量,比較與數(shù)值仿真值間誤差。通過決定系數(shù)R2、均方誤差(root mean square error, RMSE)、殘差平方和(residual sum of squares, RSS)量化誤差,擬合模型最理想誤差指標(biāo)應(yīng)為R2接近1、RMSE和RSS接近0,說明該磁體模型精度越高,擬合模型空間實(shí)際磁場分布特性越吻合。

    表1中,B是數(shù)值仿真結(jié)果,平均值為B-,不同模型計(jì)算值為B′。3個(gè)模型磁感應(yīng)場在探測平面上的MAPE值分別為21.55%、3.59%、3.94%。本文提出的混合模型是基于單個(gè)磁體模擬方法的改進(jìn),Ⅱ、Ⅲ模型的磁感應(yīng)場及分量的R2更接近1,誤差指標(biāo)RMSE、RSS較Ⅰ模型有較明顯優(yōu)化,完全滿足工程的高精度需要。作為對比,組成Ⅱ、Ⅲ模型的偶極子陣列陣元數(shù)目不同,Ⅱ模型在磁感應(yīng)場及y、z兩個(gè)方向分量上較Ⅲ模型誤差更小,因此并非一味靠增加磁偶極子數(shù)目來提高擬合精度。綜上,由橢球體與偶極子陣列組成的混合模型從精度與模型靈活度兩個(gè)方面實(shí)現(xiàn)了模型優(yōu)化:通過L-M算法確定磁矩參數(shù),較磁矩經(jīng)驗(yàn)值方法進(jìn)一步提高擬合精度,通過調(diào)整偶極子陣列陣元個(gè)數(shù),在滿足精度需要的前提下提高模型的靈活度。

    4.2" 驗(yàn)證相關(guān)參數(shù)影響

    在驗(yàn)證本文提出的混合模型提高磁場分布計(jì)算精度的基礎(chǔ)上,為避免單一尺寸模型結(jié)果的偶然性,增加不同尺寸航行器建模實(shí)驗(yàn)、不同偶極子陣元數(shù)目建模實(shí)驗(yàn),討論磁體模型中偶極子陣列陣元個(gè)數(shù)與探測距離等相關(guān)參數(shù)與建模精度的相關(guān)關(guān)系。以探測平面上磁場分布計(jì)算值與仿真值間的平均絕對百分比誤差(后文簡稱平均誤差)作為主要評(píng)價(jià)指標(biāo),分別進(jìn)行以下兩個(gè)實(shí)驗(yàn)。

    實(shí)驗(yàn) 1" 研究不同尺寸水下航行器建模影響。參考文獻(xiàn)[30]中的實(shí)驗(yàn)選取3種尺寸的航行器,按照長寬尺寸分為小、中、大3種艇型,小型艇模型參數(shù):艇長57.2 m,艇直徑8.6 m;中型艇模型參數(shù):艇長100 m,艇直徑10 m;大型艇模型參數(shù):艇長153 m,艇直徑17.3 m。為使結(jié)果具有可對比性,雖然目標(biāo)航行器尺寸改變,但保持探測高度與艇長的相對探測尺度不變,分析3種艇型磁場計(jì)算結(jié)果的探測平面高度分別為286 m、500 m、765 m。

    在對不同艇型航行器目標(biāo)的混合模型建模中,分別調(diào)整混合模型中磁偶極子陣列的陣元個(gè)數(shù)n=2,3,…,12,在對應(yīng)探測平面上由模型計(jì)算得到的磁場值與有限元仿真值間平均誤差如圖12所示。通過對不同陣元個(gè)數(shù)模型誤差進(jìn)行分析可知,混合模型的磁體模擬方法受偶極子陣列陣元數(shù)目影響,3種艇型航行器混合模型對磁場分布的計(jì)算誤差均可控制在11%以內(nèi);當(dāng)偶極子陣列陣元數(shù)目達(dá)到一定程度時(shí),建模誤差趨于穩(wěn)定,理論上可以通過調(diào)整陣元數(shù)目獲得最佳擬合模型。

    實(shí)驗(yàn) 2" 研究不同偶極子陣元數(shù)目建模影響。建立某中型艇(艇長100 m,艇直徑10 m)的磁場分布混合模型,偶極子陣列陣元數(shù)目為n=2,3,…,12,在500 m探測高度探測平面上選取相同數(shù)目探測點(diǎn)。計(jì)算不同陣元數(shù)目的混合模型磁場值,與數(shù)值仿真值間誤差與探測高度間關(guān)系如圖13所示。雖然混合模型對磁場分布的計(jì)算受探測高度,即探測距離影響較大,但在400~700 m探測范圍內(nèi)磁場分布計(jì)算值平均誤差控制在10%內(nèi),模型具有一定的精度穩(wěn)定性。保持相同探測高度,混合模型中偶極子陣列陣元數(shù)目對模型精度影響不明顯,但較多陣元數(shù)目的混合模型在探測距離較遠(yuǎn)表現(xiàn)更加穩(wěn)定。

    混合模型磁體模擬法有較高的建模實(shí)用價(jià)值。在本探測條件下建立的四元偶極子陣列混合模型,在探測點(diǎn)所在探測平面上磁感應(yīng)場計(jì)算平均誤差為3%、在探測空間最遠(yuǎn)距離處為32%。作為對比,均勻磁化單個(gè)橢球體模型在相同探測距離時(shí)磁場計(jì)算誤差分別為22%、47%,在本實(shí)驗(yàn)條件下可認(rèn)為四元偶極子陣列混合模型為最佳模型。通過多組實(shí)驗(yàn)得到以下分析結(jié)論:

    (1) 磁偶極子陣列陣元數(shù)目是影響混合模型計(jì)算精度的可變量,相關(guān)待求參數(shù)如磁矩值、陣元分布位置,同時(shí)受到陣元數(shù)目的影響。

    (2) 探測距離對混合模型空間磁場分布影響較大,分布隨著探測距離增大,空間磁場對稱特性出現(xiàn)畸變或峰值轉(zhuǎn)移。

    (3) 采用不同磁體模擬方法對磁場計(jì)算結(jié)果有較大影響。本文提出的通過L-M算法調(diào)整混合模型參數(shù)方法,比應(yīng)用經(jīng)驗(yàn)公式的單個(gè)磁體模擬方法計(jì)算更簡單、在探測距離增加時(shí)結(jié)果更有效。

    通過以上在不同尺寸航行器、不同偶極子陣列陣元數(shù)目、不同探測距離等相關(guān)參數(shù)下的實(shí)驗(yàn)理論分析,證明在實(shí)際磁異常探測技術(shù)的條件下如獲得任意探測點(diǎn)的位置及磁異常探測數(shù)據(jù),不僅可實(shí)現(xiàn)與水下航行器相對應(yīng)尺寸的混合模型空間磁場分布計(jì)算,也可通過對相關(guān)參數(shù)的調(diào)整獲得最佳擬合模型,滿足高精度應(yīng)用要求。

    5" 結(jié)" 論

    伴隨水下航行器目標(biāo)環(huán)境日益復(fù)雜,針對復(fù)雜磁源目標(biāo)的定位識(shí)別需建立更高精度磁場擬合模型。本文綜合有限元法與磁體模擬法的應(yīng)用研究,充分利用有限元法獲得復(fù)雜結(jié)構(gòu)航行器目標(biāo)的空間磁感應(yīng)強(qiáng)度分布,建立旋轉(zhuǎn)橢球體與偶極子陣列的高精度混合磁體數(shù)值模型與待求模型系數(shù)方法,討論不同航行器目標(biāo)混合陣列建模影響因素與最佳擬合模型,驗(yàn)證了該模型較好的擬合優(yōu)度與探測空間泛化性。本文討論了磁體模擬方法混合模型在全探測空間建模、模型中精確磁矩值、提升計(jì)算精度等多個(gè)方面的改進(jìn),拓展實(shí)現(xiàn)磁場數(shù)據(jù)延拓、高精度磁體模擬方法的理論應(yīng)用。

    參考文獻(xiàn)

    [1] 宋保維, 潘光, 張立川, 等. 自主水下航行器發(fā)展趨勢及關(guān)鍵技術(shù)[J]. 中國艦船研究, 2022, 17(5): 2744.

    SONG B W, PAN G, ZHANG L C, et al. Development trend and key technologies of autonomous underwater vehicles[J]. Chinese Journal of Ship Research, 2022, 17(5): 2744.

    [2] WELLER R A, LEVIN L A, GATES A R. Global observing needs in the deep ocean[J]. Frontiers in Marine Science, 2019, 6(6): 167218016.

    [3] WANG C, CUI Y, SONG X, et al. A novel underwater target detection method based on low-frequency magnetic signal generated by portable transmitter[J]. IEEE Sensors Journal, 2023, 23(8): 84598465.

    [4] SANCHEZ P, PAPAELIAS M, MARQUEZ F. Autonomous underwater vehicles: instrumentation and measurements[J]. IEEE Instrumentation and Measurement Magazine, 2020, 23(2): 105114.

    [5] 鞠建波, 郁紅波, 范趙鵬, 等. 吊放聲吶擴(kuò)展螺旋陣搜潛效能評(píng)估改進(jìn)方法[J]. 水下無人系統(tǒng)學(xué)報(bào), 2020, 28(4): 434439.

    JU J B, YU H B, FAN Z P, et al. An improved method for evaluation the submarine searching effectiveness with extension spiral array of dipping sonar[J]. Journal of Unmanned Undersea Systems, 2020, 28(4): 434439.

    [6] 胡平, 顧雪峰, 劉凱, 等. 國外航空反潛裝備現(xiàn)狀及發(fā)展趨勢[J]. 艦船電子工程, 2022, 42(9): 1012, 86.

    HU P, GU X F, LIU K, et al. Current situation and development trend of foreign aviation anti-submarine equipment[J]. Ship Electronic Engineering, 2022, 42(9): 1012, 86.

    [7] BRUNOTTE X, MEUNIER G, BONGIRAUD J. Ship magnetizations modelling by finite element method[C]∥Proc.of the IEEE Conference on Electromagnetic Field Computation, 1992.

    [8] WANG J X, REN Q. A 3-D hybrid maxwell’s equations finite-difference time-domain (ME-FDTD)/wave equation finite-element time-domain (WE-FETD) method[J]. IEEE Trans.on Antennas and Propagation, 2023, 71(6): 52125220.

    [9] CHADEBEC O, COULOMB J L, BONGIRAUD J P, et al. Recent improvements for solving inverse magnetostatic problem applied to thin shells[J]. IEEE Trans.on Magnetics, 2002, 38(2): 10051008.

    [10] IOANNIDIS G. Identification of a ship or submarine from its magnetic signature[J]. IEEE Trans.on Aerospaceamp;Electronic Systems, 1977, 13: 327329.

    [11] 曲曉慧, 楊日杰, 單志超. 潛艇磁場建模方法的分析與比較[J]. 艦船科學(xué)技術(shù), 2011, 33(3): 711.

    QU X H, YANG R J, SHAN Z C. Analysis and comparison on magnetic field modeling method of submarine[J]. Ship Science and Technology, 2011, 33(3): 711.

    [12] 林春生. 艦船磁場深度換算中的單列橢球體和磁偶極子混合模型[J]. 水雷戰(zhàn)與艦船防護(hù), 1996, 1996(3): 5458.

    LIN C S. Hybrid model of single-row ellipsoid and magnetic dipole in ship magnetic field depth conversion[J]. Mine Warfare & Ship Self-defence, 1996, 1996(3): 5458.

    [13] HOFFMANN A. Spin hall effects in metals[J]. IEEE Trans.on Magnetics, 2013, 49(10): 51725193.

    [14] 徐杰. 基于遺傳算法的潛艇三維模型及其磁場延拓方法[J]. 艦船電子工程, 2009, 29(7): 188191.

    XU J. Application of genetic algorithm in continuation of magnetic field of submarines in high areas[J]. Ship Electronic Engineering, 2009, 29(7): 188191.

    [15] 招妙妍. 基于模擬退火算法的航向?qū)崪y磁場數(shù)據(jù)建模分析[J]. 艦船電子工程, 2019, 39(9): 188190.

    ZHAO M Y. Modeling and analysis of directional field data based on simulated annealing algorithm[J]. Ship Electronic Engineering, 2019, 39(9): 188190.

    [16] NARA T, SUZUKI S, ANDO S. A closed form formula for magnetic dipole localizati-on by measurement of its magnetic field and spatial gradients[J]. IEEE Trans.on Magnetics, 2006, 42(10): 32913293.

    [17] QIN Y J, LI K Y, YAO C, et al. Magnetic anomaly detection using full magnetic gradient orthonormal basis function[J]. IEEE Sensors Journal, 2020, 20(21): 1292812940.

    [18] 張志東. 磁性材料的磁結(jié)構(gòu)、磁疇結(jié)構(gòu)和拓?fù)浯沤Y(jié)構(gòu)[J]. 物理學(xué)報(bào), 2015, 64(6): 521.

    ZHANG Z D. Magnetic structures, magnetic domains and topological magnetic textures of magnetic materials[J]. Acta Physica Sinica, 2015, 64(6): 521.

    [19] JILES D, KIARIE W. An integrated model of magnetic hysteresis, the magnetomechanical effect, and the barkhausen effect[J]. IEEE Trans.on Magnetics, 2023, 57(2): 111.

    [20] CHUNG H J, YANG C S, JUNG W J. A magnetic field separation technique for a scaled model ship through an earth’s magnetic field simulator[J]. Journal of Magnetics, 2015, 20(1): 6268.

    [21] HALGEDAH F. The dependence of magnetic domain structure upon magnetization state with emphasis upon nucleation as a mechanism for pseudo-single-domain behavior[J]. Journal of Geophysical Research, 1983: 65056522.

    [22] XIN Z, XU Q, LI Q, et al. An orthonormalized basis function based narrowband filtering algorithm for magnetic anomaly detection[C]∥Proc.of the International Congress on Image amp; Signal Processing, 2016.

    [23] 孫曉永, 張琦, 潘孟春, 等. 基于COMSOL Multiphysics的目標(biāo)磁特性仿真分析[J]. 中國測試, 2017, 43(1): 122126.

    SUN X Y, ZHANG Q, PAN M C, et al. Simulation analysis of the magnetic properties based on COMSOL multiphysics[J]. China Measurementamp;Test, 2017, 43(1): 122126.

    [24] RUBIO J D J. Stability analysis of the modified levenberg-marquardt algorithm for the artificial neural network training[J]. IEEE Trans.on Neural Networks and Learning Systems, 2020, 32(8): 35103524.

    [25] ANDREI N. Steepest descent methods[J]. Springer Optimization and its Applications, 2022, 195(1): 81107.

    [26] HANSEN E R, GREENBERG R I. An interval Newton method[J]. Applied Mathematics amp; Computation, 1983, 12(2/3): 8998.

    [27] VASIN V V. Solving nonlinear inverse problems based on the regularized modified Gauss-Newton method[J]. Doklady Mathe-matics, 2022, 105(3): 175177.

    [28] GE J, WANG S Q, DONG H B, et al. Real-time detection of moving magnetic target using distributed scalar sensor based on hybrid algorithm of particle swarm optimization and Gauss-Newton method[J]. IEEE Sensors Journal, 2020, 20(18): 1071710723.

    [29] YANG J, ZOU J. Parameter estimation of a horizontally multilayered soil with a fast evaluation of the apparent resistivity and its derivatives[J]. IEEE Access, 2020, 8: 5265252662.

    [30] 戴忠華, 張曉兵, 周穗華. 艦船磁場磁偶極子陣列模型的適用性研究[J]. 武漢理工大學(xué)學(xué)報(bào), 2020, 42(4): 99104.

    DAI Z H, ZHANG X B, ZHOU S H. Applicability of magnetic dipole array model for ship magnetic field[J]. Journal of Wuhan University of Technology, 2020, 42(4): 99104.

    作者簡介

    趙高陽(1997—),女,碩士研究生,主要研究方向?yàn)榇女惓L綔y。

    劉" 勇(1968—),男,研究員,碩士,主要研究方向?yàn)榈皖l電磁通信。

    朱平杰(1978—),男,高級(jí)工程師,主要研究方向?yàn)榈皖l電磁天線。

    向" 冰(1991—),女,高級(jí)工程師,碩士,主要研究方向?yàn)榈皖l電磁信號(hào)處理。

    周洪娟(1980—),女,副教授,博士,主要研究方向?yàn)楹Q蟠盘綔y、電波傳播。

    猜你喜歡
    磁場模型
    一半模型
    西安的“磁場”
    為什么地球有磁場呢
    重要模型『一線三等角』
    文脈清江浦 非遺“磁場圈”
    重尾非線性自回歸模型自加權(quán)M-估計(jì)的漸近分布
    《磁場》易錯(cuò)易混知識(shí)剖析
    磁場的性質(zhì)和描述檢測題
    3D打印中的模型分割與打包
    2016年春季性感磁場
    Coco薇(2016年1期)2016-01-11 16:53:24
    在线观看免费高清a一片| 真人一进一出gif抽搐免费| 国产国语露脸激情在线看| 性少妇av在线| 脱女人内裤的视频| 夜夜看夜夜爽夜夜摸 | 这个男人来自地球电影免费观看| 美女扒开内裤让男人捅视频| 在线观看舔阴道视频| 亚洲成人久久性| 日本 av在线| 久久人妻福利社区极品人妻图片| 国产区一区二久久| aaaaa片日本免费| 日本精品一区二区三区蜜桃| www日本在线高清视频| 国产高清videossex| 精品久久久久久成人av| 国产精品永久免费网站| 精品久久久久久久久久免费视频 | 又黄又粗又硬又大视频| 午夜影院日韩av| 久久香蕉国产精品| 午夜福利在线观看吧| 久久精品国产亚洲av高清一级| av电影中文网址| 男男h啪啪无遮挡| 人妻丰满熟妇av一区二区三区| 久9热在线精品视频| 电影成人av| 成人国语在线视频| 亚洲va日本ⅴa欧美va伊人久久| 亚洲专区中文字幕在线| 淫秽高清视频在线观看| 两人在一起打扑克的视频| 丁香六月欧美| 国产无遮挡羞羞视频在线观看| 久热这里只有精品99| 久久久精品欧美日韩精品| 久久午夜亚洲精品久久| av国产精品久久久久影院| 午夜激情av网站| 悠悠久久av| 欧美+亚洲+日韩+国产| 老司机靠b影院| 国产免费现黄频在线看| 超碰97精品在线观看| 亚洲人成网站在线播放欧美日韩| 99久久综合精品五月天人人| videosex国产| 美女高潮到喷水免费观看| 国产午夜精品久久久久久| 亚洲中文av在线| 欧美日韩亚洲综合一区二区三区_| 大型av网站在线播放| 后天国语完整版免费观看| 欧美成人免费av一区二区三区| 黑人操中国人逼视频| 亚洲专区中文字幕在线| av在线天堂中文字幕 | 亚洲国产欧美网| 亚洲 国产 在线| 亚洲精品国产精品久久久不卡| 制服诱惑二区| 精品一区二区三区av网在线观看| 午夜精品国产一区二区电影| 成人亚洲精品一区在线观看| 国产97色在线日韩免费| 好看av亚洲va欧美ⅴa在| 看黄色毛片网站| 天天添夜夜摸| 久久久水蜜桃国产精品网| 男女之事视频高清在线观看| 嫁个100分男人电影在线观看| 国产国语露脸激情在线看| 久久伊人香网站| 久久 成人 亚洲| 9热在线视频观看99| 成人特级黄色片久久久久久久| 国产精品1区2区在线观看.| 久久久久国产精品人妻aⅴ院| 成年版毛片免费区| 又黄又爽又免费观看的视频| 亚洲九九香蕉| 大陆偷拍与自拍| 欧美日韩乱码在线| 国产成人精品久久二区二区91| 精品久久久久久久毛片微露脸| 99国产精品99久久久久| 国产精品 欧美亚洲| 欧美大码av| 成人永久免费在线观看视频| 香蕉久久夜色| 大香蕉久久成人网| 99在线视频只有这里精品首页| 制服诱惑二区| 在线观看免费午夜福利视频| 69精品国产乱码久久久| 涩涩av久久男人的天堂| 亚洲少妇的诱惑av| 妹子高潮喷水视频| 亚洲精品在线美女| 曰老女人黄片| 又大又爽又粗| 久久人妻熟女aⅴ| 中国美女看黄片| 亚洲成a人片在线一区二区| 搡老乐熟女国产| 日韩大尺度精品在线看网址 | 久久国产精品影院| 国产成人免费无遮挡视频| 国产av精品麻豆| 亚洲九九香蕉| 久久久久国内视频| a级毛片在线看网站| 国产欧美日韩综合在线一区二区| 亚洲第一欧美日韩一区二区三区| 黑人欧美特级aaaaaa片| 色老头精品视频在线观看| 国产男靠女视频免费网站| 俄罗斯特黄特色一大片| 国产精品 欧美亚洲| 午夜两性在线视频| 天天添夜夜摸| 欧美另类亚洲清纯唯美| 999精品在线视频| av有码第一页| 欧美成人免费av一区二区三区| 97人妻天天添夜夜摸| 亚洲视频免费观看视频| 后天国语完整版免费观看| 国产成人免费无遮挡视频| 男女午夜视频在线观看| 最近最新免费中文字幕在线| 日日爽夜夜爽网站| 午夜福利欧美成人| 99riav亚洲国产免费| 免费观看精品视频网站| 久久国产乱子伦精品免费另类| 久久久精品欧美日韩精品| 亚洲黑人精品在线| 亚洲少妇的诱惑av| 久久精品aⅴ一区二区三区四区| 国产又爽黄色视频| 黄色视频不卡| 在线观看午夜福利视频| 在线观看免费午夜福利视频| 自拍欧美九色日韩亚洲蝌蚪91| 一区二区三区国产精品乱码| 一进一出抽搐gif免费好疼 | 90打野战视频偷拍视频| 18美女黄网站色大片免费观看| 国产精品亚洲av一区麻豆| 中文字幕高清在线视频| 人人澡人人妻人| 国产亚洲欧美精品永久| 91成年电影在线观看| 久久久久久久久久久久大奶| 亚洲人成电影观看| 极品人妻少妇av视频| 两个人看的免费小视频| 久久精品影院6| 午夜福利免费观看在线| 91国产中文字幕| 日韩一卡2卡3卡4卡2021年| 97超级碰碰碰精品色视频在线观看| 精品人妻在线不人妻| 日本 av在线| 国产成人精品久久二区二区免费| 国产精品电影一区二区三区| 亚洲一区二区三区不卡视频| 亚洲在线自拍视频| 欧美乱码精品一区二区三区| 操美女的视频在线观看| 淫秽高清视频在线观看| 成人国产一区最新在线观看| 久久狼人影院| a级毛片在线看网站| 91麻豆av在线| 日韩大尺度精品在线看网址 | 多毛熟女@视频| 国产1区2区3区精品| 97人妻天天添夜夜摸| 精品久久久久久,| 国产精品九九99| 最近最新免费中文字幕在线| 中出人妻视频一区二区| 日日爽夜夜爽网站| 欧美不卡视频在线免费观看 | 亚洲国产精品合色在线| 国产国语露脸激情在线看| 欧美成狂野欧美在线观看| 国产欧美日韩一区二区精品| 色综合欧美亚洲国产小说| 人人妻人人添人人爽欧美一区卜| 午夜福利影视在线免费观看| 变态另类成人亚洲欧美熟女 | 麻豆成人av在线观看| 黄色成人免费大全| 精品午夜福利视频在线观看一区| 黄片播放在线免费| 最好的美女福利视频网| 黄片播放在线免费| 97超级碰碰碰精品色视频在线观看| 一区二区三区国产精品乱码| 国产精华一区二区三区| 国产三级黄色录像| 一级黄色大片毛片| 国产精品免费视频内射| 99久久久亚洲精品蜜臀av| 久久久久国产精品人妻aⅴ院| 色老头精品视频在线观看| 国产乱人伦免费视频| 久久国产精品影院| 亚洲人成网站在线播放欧美日韩| 少妇 在线观看| 日韩免费高清中文字幕av| 亚洲av片天天在线观看| 欧美精品啪啪一区二区三区| 欧美日韩国产mv在线观看视频| 丁香六月欧美| 欧美日韩国产mv在线观看视频| 日韩国内少妇激情av| www日本在线高清视频| 亚洲专区字幕在线| 无限看片的www在线观看| 亚洲av五月六月丁香网| 一边摸一边做爽爽视频免费| 一级a爱视频在线免费观看| 中文亚洲av片在线观看爽| 亚洲 欧美 日韩 在线 免费| 亚洲精品在线美女| 12—13女人毛片做爰片一| 老熟妇仑乱视频hdxx| 18禁观看日本| 乱人伦中国视频| 99国产极品粉嫩在线观看| 欧美在线黄色| 亚洲视频免费观看视频| av有码第一页| 欧美日韩亚洲综合一区二区三区_| 日日干狠狠操夜夜爽| 免费观看人在逋| 黄频高清免费视频| 一夜夜www| 国产精品 国内视频| 久久九九热精品免费| tocl精华| 精品熟女少妇八av免费久了| 18禁黄网站禁片午夜丰满| 亚洲欧美激情综合另类| 精品久久蜜臀av无| 欧美另类亚洲清纯唯美| 免费av中文字幕在线| 亚洲在线自拍视频| 成人精品一区二区免费| 香蕉久久夜色| 国产亚洲精品一区二区www| 又黄又爽又免费观看的视频| 成人国产一区最新在线观看| 欧美av亚洲av综合av国产av| 男女下面进入的视频免费午夜 | 午夜免费成人在线视频| 欧美av亚洲av综合av国产av| 国产精品久久久人人做人人爽| 免费高清在线观看日韩| 欧美成狂野欧美在线观看| 如日韩欧美国产精品一区二区三区| 久久香蕉激情| 99热国产这里只有精品6| 中文字幕色久视频| 国产精品永久免费网站| 夜夜躁狠狠躁天天躁| 国产日韩一区二区三区精品不卡| 黑人欧美特级aaaaaa片| 日日爽夜夜爽网站| 成年版毛片免费区| 老司机靠b影院| 最好的美女福利视频网| 黄色女人牲交| 欧美激情极品国产一区二区三区| 日韩欧美免费精品| 久久精品亚洲av国产电影网| 国产精品秋霞免费鲁丝片| aaaaa片日本免费| 欧洲精品卡2卡3卡4卡5卡区| 免费一级毛片在线播放高清视频 | 法律面前人人平等表现在哪些方面| 我的亚洲天堂| 在线播放国产精品三级| 精品熟女少妇八av免费久了| 少妇的丰满在线观看| 黄色a级毛片大全视频| 超碰成人久久| 夜夜看夜夜爽夜夜摸 | 国产精品亚洲av一区麻豆| 国产激情久久老熟女| 亚洲情色 制服丝袜| 国产亚洲精品一区二区www| 欧美激情久久久久久爽电影 | 怎么达到女性高潮| 亚洲一区二区三区不卡视频| 亚洲五月色婷婷综合| 日韩免费av在线播放| 精品人妻在线不人妻| 国产成年人精品一区二区 | 国产av一区在线观看免费| 国产激情欧美一区二区| 亚洲专区字幕在线| 久久香蕉精品热| 亚洲国产精品999在线| 宅男免费午夜| 大型av网站在线播放| 操出白浆在线播放| 多毛熟女@视频| 搡老熟女国产l中国老女人| 亚洲人成网站在线播放欧美日韩| 久久欧美精品欧美久久欧美| 久久亚洲精品不卡| 成人免费观看视频高清| 交换朋友夫妻互换小说| 免费在线观看黄色视频的| svipshipincom国产片| 男男h啪啪无遮挡| 午夜两性在线视频| 国产三级在线视频| 欧美老熟妇乱子伦牲交| 天堂中文最新版在线下载| 在线观看免费午夜福利视频| 免费av中文字幕在线| 国产熟女xx| 又紧又爽又黄一区二区| 精品久久久久久,| 国产精品一区二区三区四区久久 | 一区在线观看完整版| 高清av免费在线| 色播在线永久视频| 色哟哟哟哟哟哟| 欧美黑人精品巨大| 日韩三级视频一区二区三区| 90打野战视频偷拍视频| 国产色视频综合| 国产精品久久电影中文字幕| 男男h啪啪无遮挡| a级毛片黄视频| 欧美不卡视频在线免费观看 | 日本 av在线| 精品国产亚洲在线| 成人av一区二区三区在线看| 久久国产乱子伦精品免费另类| 亚洲视频免费观看视频| 天堂中文最新版在线下载| 精品国产一区二区久久| 最新美女视频免费是黄的| 亚洲av五月六月丁香网| 日本欧美视频一区| 成人亚洲精品一区在线观看| 在线观看免费日韩欧美大片| 国产亚洲欧美98| 久久久久久久久中文| 妹子高潮喷水视频| 精品一品国产午夜福利视频| 叶爱在线成人免费视频播放| 99精品久久久久人妻精品| 国产高清国产精品国产三级| 国产精品国产高清国产av| 男女午夜视频在线观看| 成人免费观看视频高清| 午夜精品久久久久久毛片777| 国产av一区在线观看免费| 国产精品国产高清国产av| 亚洲伊人色综图| 19禁男女啪啪无遮挡网站| 日韩高清综合在线| 成年女人毛片免费观看观看9| 亚洲中文日韩欧美视频| 亚洲精品国产区一区二| 80岁老熟妇乱子伦牲交| 日日干狠狠操夜夜爽| 免费不卡黄色视频| 免费一级毛片在线播放高清视频 | 韩国av一区二区三区四区| 97碰自拍视频| 五月开心婷婷网| 国产极品粉嫩免费观看在线| 淫妇啪啪啪对白视频| 久久国产精品影院| 亚洲激情在线av| 交换朋友夫妻互换小说| 成熟少妇高潮喷水视频| 成年版毛片免费区| 日韩欧美一区二区三区在线观看| 一级,二级,三级黄色视频| 天堂影院成人在线观看| 侵犯人妻中文字幕一二三四区| 国产av一区二区精品久久| 中文字幕人妻丝袜制服| 欧美一级毛片孕妇| 热99国产精品久久久久久7| 亚洲 国产 在线| 一个人免费在线观看的高清视频| 日韩一卡2卡3卡4卡2021年| 午夜福利影视在线免费观看| 麻豆一二三区av精品| 亚洲 欧美 日韩 在线 免费| 天堂影院成人在线观看| 国产精品乱码一区二三区的特点 | 亚洲人成网站在线播放欧美日韩| 国产av一区二区精品久久| 视频在线观看一区二区三区| 国产精品永久免费网站| 精品国产超薄肉色丝袜足j| 亚洲成人国产一区在线观看| 成人黄色视频免费在线看| 免费在线观看日本一区| 色哟哟哟哟哟哟| 亚洲精品美女久久av网站| 欧美色视频一区免费| 最近最新中文字幕大全电影3 | av在线播放免费不卡| 99精国产麻豆久久婷婷| 男人操女人黄网站| 欧美大码av| 九色亚洲精品在线播放| 国产一卡二卡三卡精品| 久久久久九九精品影院| 两个人看的免费小视频| 亚洲欧美日韩另类电影网站| 亚洲精品一二三| 国产精品偷伦视频观看了| 国产成人av激情在线播放| 又紧又爽又黄一区二区| 久久天堂一区二区三区四区| 丝袜人妻中文字幕| 国产人伦9x9x在线观看| 在线观看66精品国产| 免费在线观看完整版高清| 国产色视频综合| 一级毛片高清免费大全| 国产精品影院久久| 少妇裸体淫交视频免费看高清 | 黄色女人牲交| www.www免费av| 亚洲精品av麻豆狂野| 午夜影院日韩av| 亚洲精华国产精华精| 91大片在线观看| 免费在线观看影片大全网站| 琪琪午夜伦伦电影理论片6080| 首页视频小说图片口味搜索| 国产精品99久久99久久久不卡| a级片在线免费高清观看视频| 18禁裸乳无遮挡免费网站照片 | 国产免费男女视频| 少妇被粗大的猛进出69影院| 国产精品自产拍在线观看55亚洲| 亚洲aⅴ乱码一区二区在线播放 | 最近最新中文字幕大全免费视频| 99久久综合精品五月天人人| 国产精品1区2区在线观看.| av欧美777| 亚洲五月天丁香| 窝窝影院91人妻| 80岁老熟妇乱子伦牲交| 人人妻人人澡人人看| 欧美乱色亚洲激情| 日韩精品免费视频一区二区三区| 欧美日韩乱码在线| 免费久久久久久久精品成人欧美视频| 9191精品国产免费久久| 亚洲 国产 在线| 国产伦人伦偷精品视频| 欧美中文综合在线视频| 999精品在线视频| 午夜福利免费观看在线| 国产乱人伦免费视频| 国产精品国产av在线观看| 欧美乱妇无乱码| 久久久久久久久中文| 一个人观看的视频www高清免费观看 | a级毛片黄视频| 久久精品国产亚洲av香蕉五月| 超色免费av| 久久久久久久久中文| 精品一品国产午夜福利视频| 人人妻人人澡人人看| 亚洲精品成人av观看孕妇| 国产有黄有色有爽视频| 一级片'在线观看视频| 午夜两性在线视频| 午夜福利在线免费观看网站| 国产成人啪精品午夜网站| 日韩中文字幕欧美一区二区| 琪琪午夜伦伦电影理论片6080| 成年女人毛片免费观看观看9| 在线观看www视频免费| 亚洲人成77777在线视频| 午夜成年电影在线免费观看| 亚洲专区中文字幕在线| 91av网站免费观看| 亚洲在线自拍视频| 涩涩av久久男人的天堂| 十八禁人妻一区二区| av电影中文网址| 欧美黄色片欧美黄色片| 久久久久久亚洲精品国产蜜桃av| 老汉色∧v一级毛片| 欧美日韩亚洲综合一区二区三区_| 成人影院久久| 最好的美女福利视频网| 亚洲人成电影观看| 精品国产国语对白av| 少妇被粗大的猛进出69影院| 老司机午夜十八禁免费视频| 亚洲专区字幕在线| 久久天堂一区二区三区四区| 极品人妻少妇av视频| 在线观看www视频免费| 免费女性裸体啪啪无遮挡网站| 午夜免费鲁丝| 一级a爱片免费观看的视频| 别揉我奶头~嗯~啊~动态视频| 国产精品乱码一区二三区的特点 | 欧美+亚洲+日韩+国产| 电影成人av| 国产亚洲精品综合一区在线观看 | 国产高清激情床上av| 国产麻豆69| 久久精品91蜜桃| а√天堂www在线а√下载| xxxhd国产人妻xxx| 两个人免费观看高清视频| 夜夜夜夜夜久久久久| 悠悠久久av| www.www免费av| 日韩免费高清中文字幕av| 久久久水蜜桃国产精品网| 美女扒开内裤让男人捅视频| 热re99久久精品国产66热6| 久久99一区二区三区| 在线国产一区二区在线| 午夜免费观看网址| 成人三级做爰电影| 精品人妻1区二区| 国产精品偷伦视频观看了| 欧美不卡视频在线免费观看 | 午夜精品国产一区二区电影| 电影成人av| 日韩国内少妇激情av| 欧美日韩中文字幕国产精品一区二区三区 | 国产无遮挡羞羞视频在线观看| 免费人成视频x8x8入口观看| 精品国产国语对白av| 国产黄a三级三级三级人| 欧美大码av| 在线观看免费视频网站a站| 色婷婷av一区二区三区视频| 动漫黄色视频在线观看| 国产精品 国内视频| 香蕉久久夜色| 19禁男女啪啪无遮挡网站| 黄色a级毛片大全视频| 欧美午夜高清在线| av网站在线播放免费| 涩涩av久久男人的天堂| 久久九九热精品免费| 国产精品99久久99久久久不卡| 国产精品爽爽va在线观看网站 | 国产熟女xx| 女生性感内裤真人,穿戴方法视频| av天堂久久9| 国产真人三级小视频在线观看| 久久久精品国产亚洲av高清涩受| ponron亚洲| 69精品国产乱码久久久| 久久国产精品影院| 日韩人妻精品一区2区三区| 午夜91福利影院| 99在线人妻在线中文字幕| 欧美丝袜亚洲另类 | 最近最新免费中文字幕在线| 一本大道久久a久久精品| 视频区欧美日本亚洲| 欧美久久黑人一区二区| 一进一出好大好爽视频| 久久久国产一区二区| av中文乱码字幕在线| 午夜精品久久久久久毛片777| 在线观看免费视频日本深夜| 黄色视频,在线免费观看| 伦理电影免费视频| 欧美一区二区精品小视频在线| 丝袜美腿诱惑在线| 国产片内射在线| 国产精品综合久久久久久久免费 | 久久久久国产一级毛片高清牌| 日韩中文字幕欧美一区二区| 国产精品久久视频播放| www日本在线高清视频| 亚洲精品中文字幕一二三四区| 欧美 亚洲 国产 日韩一| 极品教师在线免费播放| 在线国产一区二区在线| 亚洲精品国产精品久久久不卡| 精品国内亚洲2022精品成人| 99国产精品免费福利视频| 黑人欧美特级aaaaaa片| 男女床上黄色一级片免费看| 久久人人爽av亚洲精品天堂| 久久中文看片网| 国产日韩一区二区三区精品不卡| 国产一卡二卡三卡精品| 中文字幕最新亚洲高清| 亚洲视频免费观看视频| 熟女少妇亚洲综合色aaa.|