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

    二階傳輸條件撕裂互連法在電大輻射中的應用

    2015-11-17 21:16:34孔才智郝毫保
    現(xiàn)代電子技術 2015年16期
    關鍵詞:有限元法電磁輻射

    孔才智+郝毫保

    摘 要: 針對于大尺寸電磁輻射問題,將撕裂互連法應用于三維電大尺寸輻射問題的計算仿真。該算法是區(qū)域分解方法中的一種,求解區(qū)域劃分成互不重疊的子區(qū)域,各子區(qū)域之間通過拉格朗日乘子將交界面的連續(xù)邊界條件或者傳輸條件耦合。鑒于原來的傳輸條件存在不收斂或者收斂較慢情況,提出同時能在TE/TM凋落模式快速收斂的二階傳輸邊界條件作為子區(qū)域之間交界面的邊界條件,并且引入虛擬激勵流為交換信息。數(shù)值結果表明,該改進算法有效地提高了撕裂互連法在迭代求解中的收斂性。在求解三維電大輻射問題時,該算法仿真結果與有限元算法結果一致,表明撕裂互連算法是計算大尺寸電磁輻射問題的一種有效方法。

    關鍵詞: 撕裂互連法; 有限元法; 二階傳輸邊界條件; 電磁輻射

    中圖分類號: TN911?34; O441 文獻標識碼: A 文章編號: 1004?373X(2015)16?0004?05

    Application of tearing and interconnecting method in large?scale electromagnetic radiation under second?order transmission conditions problem

    KONG Caizhi, HAO Haobao

    (South China Normal University, Guangzhou 510006, China)

    Abstract: The tearing and interconnecting method is applied to the simulation of three?dimensional (3?D) large?scale electromagnetic radiation to solve the problems of large?scale electromagnetic radiation. As one of the domain decomposition methods, this method reduces the computational complexity by dividing the domain under solution into several non?overlapping subdomains. The continuous boundary conditions or transmission conditions of the sub?domain interfaces are coupled by using the Lagrange multipliers and an iterative process is used to realize the solution of the original global problem. For further increasing the convergence of this method, a second?order transmission condition that can transmit both TE and TM evanescent modes is developed A virtual excitation source is introduced as the exchange information. Numerical simulation results reveal that the modification has improved effectively convergence performance in the iteration solution of the tearing and interconnecting method. Regarding to the electromagnetic radiation problems, the results are well matched with those obtained by finite element method.

    Keywords: tearing and interconnecting method; finite element method; second?order transmission condition; electromagnetic radiation

    0 引 言

    傳統(tǒng)的有限元方法在處理電大尺寸輻射問題和復雜細節(jié)結構目標上出現(xiàn)未知數(shù)的數(shù)量很多,導致全局線性系統(tǒng)矩陣規(guī)模龐大,矩陣形態(tài)較差,迭代求解往往很慢甚至不收斂。區(qū)域分解法提出在并行計算上的優(yōu)勢,是有限元方法解決這些問題的重要方法之一。在眾多區(qū)域分解算法當中,撕裂互連法被證明是非常有前景和高效的算法之一[1?3]。該算法將原求解區(qū)域劃分成若干個不重疊的子區(qū)域,在子區(qū)域的交界面上使用拉格朗日乘子隱式強加場的連續(xù)性條件把子區(qū)域連接起來,通過消去法將全局未知量的矩陣變成子區(qū)域交界面未知量的系統(tǒng)矩陣,一旦獲得交界面上的未知量,代回到相應的子區(qū)域,各子區(qū)域就可以單獨求解。撕裂互連算法不僅具有區(qū)域分解法的“分而治之”的策略優(yōu)點,而且減少全局未知量為交界面的未知量,使其系統(tǒng)矩陣應用Krylov子空間迭代求解[4?5],避免了求解大矩陣的復雜問題。為了提高算法在高頻電磁中的迭代收斂性,實現(xiàn)子區(qū)域交界面的一階Robin類型傳輸條件(FOTC)替代連續(xù)性條件耦合,并將該算法應用到三維電大尺寸電磁問題上,模型用一階吸收邊界條件或者完全匹配層(PML)截斷求解區(qū)域[6?8]。由于一階傳輸條件(FOTC)只能用于傳輸模式,為了能夠傳輸TE的凋落模式和提高迭代算法的收斂性,在子區(qū)域之間的交界面引入二階傳輸條件(SOTC?TE)求解大規(guī)模多區(qū)域的電磁模型[9?10]。具體的計算實例表明,二階傳輸條件具有更高的收斂效果。

    為了達到提高算法精度和分析復雜結構的目標,網(wǎng)格剖分變密,導致分界面的單元中存在凋落模式,凋落模式會嚴重影響到傳輸條件的收斂模式,上述提到的二階傳輸條件(SOTC?TE)僅能在TE模式下快速收斂。本文將推導出適用于撕裂互連法新型的SOTC?FULL二階傳輸條件,使得該傳輸條件具備更明顯的物理意義,該傳輸條件能同時在TE和TM凋落模式下快速收斂。將改進后的算法應用于求解三維電磁輻射問題中,采用一階吸收邊界條件截斷求解區(qū)域。

    1 撕裂互連原理與傳輸條件

    一個有截斷邊界的有源電磁問題可以歸結為在[Ω]空間上求解邊值問題的矢量波動方程。不失一般性,將求解自由空間劃分為互不重疊的子空間。每個子空間內建立有限元矢量波動方程以及限定的邊界條件:

    在[Ωs]子空間內滿足:

    [?×μ-1r?×E-k2oεrE=-jkozoJimp] (1)

    在PEC邊界滿足的邊界條件:

    [n×E=0] (2)

    在PMC邊界滿足的邊界條件:

    [n×?×E=0] (3)

    一階吸收邊界條件(ABC):

    [n×(?×E)+jkin×(n×E)=Uinc] (4)

    其中:[Ωs]表示第[s]子空間;[ko]和[zo]分別為自由空間的波數(shù)和波阻抗;[Jimp]為內部強加源;PEC為理想導體邊界條件;PMC為理想磁場邊界條件;[Uinc=][n×(?×Einc)+jkin×(n×Einc)]為入射波,[ki=koμrεr]。

    當一個原求解區(qū)域劃分為若干個子區(qū)域,只要保證子區(qū)域之間場的連續(xù)性,則由子區(qū)域與子區(qū)域之間的傳輸條件所構成的問題與原問題等效,這是區(qū)域分解法的基本思想。撕裂互連法是將拉格朗日乘子引入到傳輸條件中,最初的算法是用Dirichlet邊界連續(xù)性條件表示傳輸條件,由于在該條件下存在無法收斂或者收斂很慢的情況,因此引用Robin一階傳輸條件;進一步應用發(fā)現(xiàn),當網(wǎng)格劃分變密,頻率變低時[11?12],分界面單元中將存在凋落模式,該模式影響傳輸條件的收斂效果。為了改善收斂效果,本文引入非共形區(qū)域分解法(Non?conformal Domain Decomposition Method)[13?15]中能同時在TE和TM凋落模式下快速收斂的二階傳輸條件(SOTC?FULL),推導出使該方法適用于撕裂互連算法的形式,表示如下:

    [ni×(μ-1ri?×Ei)+αini×(ni×Ei)-βi?× [ni(?×Ei)n]+γi?τ?τ·[ni×(μ-1ri?×Ei)]=Λ] (5)

    式中:[i]表示子區(qū)域之間的界面;[αi=jkoμsqrεsqr],[μsqr]和[εsqr]分別表示子區(qū)域界面兩邊的相對介電常數(shù)和磁導率的均值;[βi]和[γi]由最小網(wǎng)格和基函數(shù)的階數(shù)決定,[βi=-j(k+ktez)],[γi=1(k2+kktmz)],在這里取[kmax,ter=πhmin],[kmax,tmr=π/(2×hmin)],[hmin]是子區(qū)域界面的最小網(wǎng)格尺寸,若[βi=0],[γi=0],則式(5)為一階傳輸條件FOTC,若[γi=0],式(5)變?yōu)镾OTC?TE傳輸條件,[Λ]為子區(qū)域界面引入的拉格朗日乘子未知量。

    基于上述邊界條件的基礎,本文推導出一種新型的二階傳輸條件,對于式(5)中的第4項[γi?τ?τ·[ni×(μ-1ri?×Ei)]],因為[ni×(μ-1ri?×Ei)]只有切向分量,所以:

    [?τ·[n×(μ-1r?×E)]=(?-?n)·[n×(μ-1r?×E)] =[?·[n×(μr-1?×E)]]]

    根據(jù)矢量恒等式和在自由空間的無源二階矢量波動方程得到:

    [?·[n×(μ-1?×E)]=n·[?×(μ-1?×E)]=k2oεrn·E]

    將式(5)中的傳輸條件展開分別表示為交界面2個不同子區(qū)域的傳輸條件,把分界面上的傳輸條件看成是一種激勵,本文中以[gi]來表示,這一激勵類似于面電流。因此式(5)可以分解為以下2種形式:

    [ni1×(μ-1ri1?×Ei1)+αi1ni1×(ni1×Ei1)-βi1?× [ni1(?×Ei1)n]+γi1?τn·E=gi] (6)

    [-ni2×(μ-1ri2?×Ei2)+αi2ni2×(ni2×Ei2)-βi2?× [ni2(?×Ei2)n]-γi2?τn·E=gi] (7)

    式中:[i1]表示[Si∈Ω1];[i2]表示[Si∈Ω2]。

    對于子區(qū)域[Ωs]聯(lián)立方程式(1)~式(4),式(6),式(7),將子區(qū)域矢量場用矢量棱邊元展開,并運用伽遼金方法,獲得子區(qū)域的矩陣方程:

    [Ks{Εs}={fs}-i∈neibourΩssiNs·gids] (8)

    式中:s表示子區(qū)域編號;[Εs]為待求的電場未知量;子區(qū)域有限元矩陣[Ks]和激勵向量[fs]的具體形式如下:

    [Ks=Ω[μ-1r(?×Ns)·(?×Ns)-k2oεrNs·Ns]dv+ jkoSABC(n×Ns)·(n×Ns)ds+αs(n×Ns)·(n×Ns)ds+ βs(?×Ns)n(?×Ns)nds+γs(n·Ns)(?τ·Ns)ds] (9)

    [{fs}=-jkozoΩNs·Jimpdv-SABCNs·Uincds] (10)

    每個子區(qū)域中的單元邊未知量按照位置可以分為3種類型[E=[Ev,Ei,Ec]],分別表示內部棱邊,交界面和拐角處未知量,其中[Ev,Ei]可以歸為一類[Er]。因此子區(qū)域的矩陣有限元系統(tǒng)可以表示為:

    [Ks=KsrrKsrcKscrKscc, {Εs}=ΕsrΕsc] (11)

    交界面[Si]上定義與交界面相關的矢量[ρi]以及布爾矩陣[Bi]。假定某一棱邊在子區(qū)域[Ωi]和子區(qū)域[Ωj]的交界面[Si]上的編號[p],在子區(qū)域內[Ωi]的編號為[q],在子區(qū)域[Ωj]的編號為[t],則[BTi(p,q)=1],[BTj(p,t)=1],因此:

    [BTiρi=-BTjρi=siNs·gids] (12)

    將所有的拐角(多個子區(qū)域重疊的邊)處的棱邊進行編號,并定義布爾矩陣[Bsc]為棱邊對應在拐角上的編號與子區(qū)域內部編號的映射關系矩陣。假定在子區(qū)域內[Ωs]內某一拐角的編號為[m],在全局拐角棱邊編號為[n],則[Bsc(m,n)=1],用[Ec]表示全局拐角未知量:

    [BscEc=Esc] (13)

    聯(lián)立式(8),式(11)~式(13),在子區(qū)域[Ωs]內矩陣形式可表示為:

    [KsrrKsrcKscrKsccΕsrΕsc=fsrfsc--Biρi+Bjρi 0] (14)

    不失一般性,將整個求解區(qū)域劃分為4個子區(qū)域,如圖1所示。

    圖1 區(qū)域分解圖

    在子區(qū)域之間交界面[S1],子區(qū)域[Ω1]和[Ω2]之間的界面,同時顯示強加Dirichlet連續(xù)性邊界條件[B1E1+B2E2=0]。整合4個子區(qū)域系統(tǒng)矩陣方程式(14),得到全局系統(tǒng)矩陣[K]如下:

    由式(16)可以看出,本文推導出一種新型的二階傳輸條件以及引入的虛擬激勵流后,獲得的系統(tǒng)矩陣與原有撕裂互連算法中的系統(tǒng)矩陣一致,很好地保證了撕裂互連算法區(qū)別于其他區(qū)域分解法最大特點[5],因此可以應用撕裂互連算法中的求解方法求解系統(tǒng)線性矩陣方程。應用高斯消去法,消去全局拐角未知量[Ec],可以獲得關于交界面上的對偶未知量[λ]的全局系統(tǒng)方程為:

    [[Krr+KrcK-1ccKTrc]λ=dr-KrcK-1ccfc] (17)

    [{Ec}=[Kcc]-1({fc}+[Krc]T{λ})] (18)

    [{Esr}=[Ksrr]-1({fsr}-[Ps]T{λ}-[Ksrc][Bsc]{Ec})] (19)

    式中[Krr],[Krc],[Kcc],[Krc],[dr],[fc]的具體表達式以及以上公式的詳細推導過程見文獻[16]。式(17)稱為全局界面系統(tǒng)方程,未知量[λ]只定義在區(qū)域之間的界面。經(jīng)過驗證,該方程求解通過Krylov子空間迭代方法如穩(wěn)態(tài)雙共軛梯度算法(BiGCSTAB)和廣義最小余量算法(GMRES)求解更有效。獲得交界面的未知量[λ]后,各個子域內部的電場值可以通過子域方程式(18),式(19)分別求出。原來三維有限元方程問題轉化為二維界面問題,交界面上的未知量通常比較少,因此求解該方程的維度和計算量比原有限元系統(tǒng)方程減少很多,所以該方法能夠求解電大電磁輻射和散射問題。

    2 數(shù)值算例

    算例1:為了驗證撕裂互連法在不同的傳輸條件FOTC,SOTC?TE,SOTC?FULL下的收斂性和算法的正確性,仿真了自由空間波傳播的算例,如圖2所示。

    圖2 自由空間波傳播求解區(qū)域

    該模型內外都是空氣,入射波[E=xe-jkz+ye-jkz],入射波沿Z軸傳播的與X軸成45°極化方向的平面波。求解空間用一階吸收邊界條件截斷外邊界,整體求解區(qū)域剖分為4個子區(qū)域,網(wǎng)格劃分密度為[λ20],入射波頻率為10 GHz,收斂精度為[10-8],采用完全相同的網(wǎng)格剖分和子區(qū)域劃分方式,全局具有2 240個未知量。從圖3迭代結果可以看出SOTC?FULL的收斂速度比FOTC,SOTC?TE傳輸條件的收斂效果好。

    圖3 自由空間波傳播模型不同傳輸條件迭代結果對比圖

    算例2:因為撕裂互連法適用于大尺寸電磁分析特別是天線有限陣列,因為在多數(shù)典型的天線陣列中,很多天線單元是相同的。對天線陣列進行仿真,只要構造和計算幾個子區(qū)域的有限元矩陣。整個求解區(qū)域劃分為9個子區(qū)域,如圖4所示。

    圖4 具有不同有限元矩陣的9個子區(qū)域圖

    由于應用外邊界如吸收邊界條件和完全匹配層,邊和角的有限元矩陣不同于內部的有限元矩陣,為此,其他的有限元矩陣可由這9個矩陣確定下來。也就是說,這9個子區(qū)域的有限元矩陣[Ksrr]一旦構成和計算出來,這些矩陣可以在全局系統(tǒng)界面方程迭代過程和最后的各個子區(qū)域求解過程中替代其他子區(qū)域矩陣。

    本文應用撕裂對接算法仿真和分析微帶饋電的貼片天線,如圖5所示。

    圖5 微帶饋電貼片天線結構圖

    16 mm12.45 mm的金屬貼片位于0.794 mm厚,相對介電常數(shù)[εr=2.2]的基底,微帶寬度為2.46 mm,特征阻抗為50 Ω,微帶與貼片的距離為2.09 mm。求解區(qū)域用吸收邊界條件截斷求解空間,激勵方法為電流激勵[J=zIδx-x0δy-y0]。

    圖6給出了撕裂互連算法求解該模型的S11參數(shù)結果圖,文章中的算法計算結果和ANSYS求解的結果吻合的非常好。

    圖6 貼片天線的S11參數(shù)算法仿真與ANSYS仿真結果對比圖

    3 結 語

    本文推導和提出適用于撕裂互連算法的二階傳輸條件(SOTC?FULL),并簡要介紹了其公式推導過程。驗證了該二階傳輸條件(SOTC?FULL)比一階傳輸條件(FOTC)有更好的收斂性,并且將其應用到撕裂互連算法求解電磁輻射問題上分析微帶饋電貼片天線。仿真結果與ANSYS軟件中的求解結果符合程度高,從而證明了撕裂互連算法是計算大尺寸電磁輻射問題的一種有效方法。

    參考文獻

    [1] LI Y, JIN J M. A vector dual?primal finite element tearing and interconnecting method for solving 3?D large?scale electromagnetic problems [J]. IEEE Transactions on Antennas and Propagation, 2006, 54(10): 3000?3009.

    [2] LI Y J, JIN J M. A new dual?primal domain decomposition approach for finite element simulation of 3?D large?scale electromagnetic problems [J]. IEEE Transactions on Antennas and Propagation, 2007, 55(10): 2803?2810.

    [3] 杜磊.時域有限元電磁計算方法的研究[D].南京:南京理工大學,2010.

    [4] ZHAO Kezhong. A domain decomposition method for solving electrically large electromagnetic problems [D]. Ohio State: The Ohio State University, 2007.

    [5] JIN J M, RILEY D J. Finite element analysis of antennas and arrays[M]. John Wiley & Sons, 2009. 336?387.

    [6] 賈會亮.頻域有限元區(qū)域分解法的研究[D].南京:南京理工大學,2009.

    [7] XUE M F, JIN J M. Nonconformal FETI?DP methods for large?scale electromagnetic simulation [J]. IEEE Transactions on Antennas and Propagation, 2012, 60(9): 4291?4305.

    [8] 高紅偉,盛新慶.三維非均勻電大目標散射的有限元撕接區(qū)域分解法計算[J].北京理工大學學報,2012,32(6):602?606.

    [9] PENG Zhen, RAWAT Vineet, LEE Jin?Fa. One way domain decomposition method with second order transmission conditions for solving electromagnetic wave problesms [J]. Journal of Computational Physics, 2010, 229: 1181?1197.

    [10] XUE M F, JIN J M. A hybrid conformal/nonconformal domain decomposition method for multi?region electromagnetic modeling[J]. IEEE Transactions on Antennas and Propagation, 2014, 62: 2009?2021.

    [11] 馬金.非均勻介質中低頻近場的分析與探測應用[D].成都:電子科技大學,2013.

    [12] YAO W, JIN J M, KREIN P T. Application of the LU recombination method to the FETI?DP method for solving low?frequency multiscale electromagnetic problems [J]. IEEE Transactions on Magnetics, 2013, 49(10): 5346?5355.

    [13] PENG Z, LEE J F. Non?conformal domain decomposition method with mixed true second order transmission condition for solving large finite antenna arrays [J]. IEEE Transactions on Antennas and Propagation, 2011, 59(5): 1638?1651.

    [14] PENG Z, LEE J F. A scalable nonoverlapping and nonconformal domain decomposition method for solving time?harmonic maxwell equations in R^3 [J]. SIAM Journal on Scientific Computing, 2012, 34(3): A1266?A1295.

    [15] MA J, NIE Z, SUN X. Efficient modeling of large?scale electromagnetic well?logging problems using an improved nonconformal FEM?DDM [J]. IEEE Transactions on Geoscience and Remote Sensing, 2014, 52(3): 1825?1833.

    [16] YAO W. Accurate, efficient, and stable domain decomposition methods for analysis of electromechanical problems [D]. Urbana: University of Illinois at Urbana?Champaign, 2014.

    猜你喜歡
    有限元法電磁輻射
    4.微波爐的電磁輻射大嗎?
    正交各向異性材料裂紋疲勞擴展的擴展有限元法研究
    論電磁輻射對環(huán)境的污染及防護措施
    活力(2019年15期)2019-09-25 07:21:58
    基于某型車載電臺的電磁輻射研究及整改
    電子測試(2018年22期)2018-12-19 05:12:12
    移動衛(wèi)星通信車電磁輻射分析與防護
    三維有限元法在口腔正畸生物力學研究中發(fā)揮的作用
    集成對稱模糊數(shù)及有限元法的切削力預測
    電磁輻射儀在確定上覆煤柱及采空區(qū)探測中的應用
    有限元法在機械設計方向中的教學實踐
    基于HCSR和CSR-OT的油船疲勞有限元法對比分析
    船海工程(2013年6期)2013-03-11 18:57:25
    国产亚洲精品久久久com| 成人二区视频| 国产69精品久久久久777片| 少妇人妻精品综合一区二区| 国产黄色视频一区二区在线观看| 老司机影院成人| 狂野欧美白嫩少妇大欣赏| 国产视频内射| 国产精品国产三级国产专区5o| 亚洲欧美日韩另类电影网站 | 在线观看一区二区三区激情| 在线观看一区二区三区激情| 在线观看免费视频网站a站| 国产欧美日韩精品一区二区| 国产淫语在线视频| 国产成人精品一,二区| 在线观看免费视频网站a站| 色综合色国产| 国产伦精品一区二区三区视频9| 亚洲精品久久午夜乱码| 久久人人爽人人片av| tube8黄色片| 亚洲最大成人中文| 精品国产三级普通话版| 夜夜骑夜夜射夜夜干| 精品亚洲乱码少妇综合久久| 伦理电影大哥的女人| 亚洲天堂av无毛| 丰满迷人的少妇在线观看| 亚洲人与动物交配视频| 少妇的逼水好多| 欧美精品一区二区大全| 有码 亚洲区| 大片免费播放器 马上看| 丰满少妇做爰视频| 人妻少妇偷人精品九色| 国产成人aa在线观看| 国产精品国产三级国产专区5o| 国产69精品久久久久777片| 久久精品国产亚洲av天美| 中文字幕久久专区| 久久久欧美国产精品| 看免费成人av毛片| 精品熟女少妇av免费看| 午夜福利在线观看免费完整高清在| 久久热精品热| 久久99蜜桃精品久久| 五月玫瑰六月丁香| 国产黄色视频一区二区在线观看| 人妻 亚洲 视频| av网站免费在线观看视频| 亚洲av福利一区| 2021少妇久久久久久久久久久| videos熟女内射| 久久久久久久久久人人人人人人| 亚洲怡红院男人天堂| 国产精品久久久久久av不卡| 国产爽快片一区二区三区| 国产黄片美女视频| av播播在线观看一区| 国产视频首页在线观看| 免费看日本二区| 亚洲在久久综合| 精品人妻视频免费看| 天堂8中文在线网| 少妇熟女欧美另类| 看免费成人av毛片| 国产精品av视频在线免费观看| 久久国产亚洲av麻豆专区| 女的被弄到高潮叫床怎么办| 成人高潮视频无遮挡免费网站| av国产免费在线观看| 亚洲色图综合在线观看| 免费黄频网站在线观看国产| 你懂的网址亚洲精品在线观看| 亚洲一级一片aⅴ在线观看| 中文字幕精品免费在线观看视频 | 高清av免费在线| 简卡轻食公司| 欧美日本视频| 日本色播在线视频| 97热精品久久久久久| 精品亚洲成国产av| 在线观看av片永久免费下载| 精品一区在线观看国产| a级毛色黄片| 99re6热这里在线精品视频| 精品人妻视频免费看| 亚洲,一卡二卡三卡| 中文字幕制服av| 岛国毛片在线播放| 亚洲精品色激情综合| 在线 av 中文字幕| 麻豆成人av视频| 久久精品久久久久久噜噜老黄| 毛片一级片免费看久久久久| 最近最新中文字幕大全电影3| av.在线天堂| 午夜激情久久久久久久| 国产乱来视频区| 国产精品99久久久久久久久| 婷婷色综合大香蕉| 久久国内精品自在自线图片| 男女无遮挡免费网站观看| 1000部很黄的大片| 超碰av人人做人人爽久久| 极品教师在线视频| .国产精品久久| 日日撸夜夜添| 免费黄色在线免费观看| 五月开心婷婷网| 久久久久久伊人网av| 日韩免费高清中文字幕av| 亚洲av不卡在线观看| 亚洲色图av天堂| 亚洲综合色惰| 一级黄片播放器| 亚洲性久久影院| 日韩伦理黄色片| videossex国产| 成人美女网站在线观看视频| 91午夜精品亚洲一区二区三区| 99九九线精品视频在线观看视频| 美女视频免费永久观看网站| 网址你懂的国产日韩在线| 日日撸夜夜添| 狂野欧美白嫩少妇大欣赏| 久久久久久久久久成人| 菩萨蛮人人尽说江南好唐韦庄| 婷婷色综合www| 亚洲丝袜综合中文字幕| av播播在线观看一区| 久久 成人 亚洲| 亚洲人成网站在线播| 在线天堂最新版资源| 国产国拍精品亚洲av在线观看| 免费黄色在线免费观看| 亚洲成人av在线免费| 五月玫瑰六月丁香| 特大巨黑吊av在线直播| 精品国产露脸久久av麻豆| 自拍欧美九色日韩亚洲蝌蚪91 | 国产成人午夜福利电影在线观看| 日本午夜av视频| 国产免费一区二区三区四区乱码| 成人综合一区亚洲| 欧美日韩精品成人综合77777| 97热精品久久久久久| 国产欧美日韩精品一区二区| 成人黄色视频免费在线看| 王馨瑶露胸无遮挡在线观看| 成人一区二区视频在线观看| 亚洲欧美日韩无卡精品| 欧美成人一区二区免费高清观看| 免费人成在线观看视频色| 一区二区三区精品91| 少妇裸体淫交视频免费看高清| 麻豆成人午夜福利视频| 国产大屁股一区二区在线视频| 女人久久www免费人成看片| 亚洲一级一片aⅴ在线观看| 少妇人妻精品综合一区二区| 色综合色国产| 日韩强制内射视频| 国产日韩欧美亚洲二区| 久久99蜜桃精品久久| 18禁在线无遮挡免费观看视频| 国产日韩欧美亚洲二区| 成年美女黄网站色视频大全免费 | 自拍欧美九色日韩亚洲蝌蚪91 | 日本黄色日本黄色录像| 亚洲av福利一区| 美女高潮的动态| 欧美最新免费一区二区三区| 亚洲欧洲日产国产| 日日啪夜夜爽| 欧美精品一区二区免费开放| 成人无遮挡网站| 小蜜桃在线观看免费完整版高清| 少妇被粗大猛烈的视频| 热99国产精品久久久久久7| 国产午夜精品一二区理论片| 久久精品国产自在天天线| 日本一二三区视频观看| 黄色配什么色好看| 精品一区在线观看国产| 亚洲国产毛片av蜜桃av| 国产一区二区三区av在线| 亚洲国产最新在线播放| 久久亚洲国产成人精品v| 婷婷色麻豆天堂久久| 国产视频首页在线观看| 成人二区视频| 欧美激情极品国产一区二区三区 | 中国美白少妇内射xxxbb| 精品一区二区免费观看| 美女高潮的动态| 亚洲精品成人av观看孕妇| 成人高潮视频无遮挡免费网站| 精品久久久久久电影网| 久久精品夜色国产| 欧美最新免费一区二区三区| 激情 狠狠 欧美| 涩涩av久久男人的天堂| 欧美日韩视频精品一区| 久久久久精品性色| 最近中文字幕2019免费版| 一级毛片电影观看| 亚洲图色成人| 欧美日韩视频精品一区| 97精品久久久久久久久久精品| 少妇高潮的动态图| 一级二级三级毛片免费看| 观看免费一级毛片| 亚洲欧洲国产日韩| 美女高潮的动态| 亚洲电影在线观看av| 舔av片在线| 97在线人人人人妻| av国产免费在线观看| av女优亚洲男人天堂| 人妻制服诱惑在线中文字幕| 亚洲内射少妇av| 亚洲人成网站在线播| 日韩一本色道免费dvd| 欧美成人a在线观看| 80岁老熟妇乱子伦牲交| 最新中文字幕久久久久| 色视频www国产| 欧美三级亚洲精品| 天天躁日日操中文字幕| 免费看不卡的av| 日日撸夜夜添| 蜜桃久久精品国产亚洲av| 日韩av不卡免费在线播放| 欧美xxxx黑人xx丫x性爽| 日日摸夜夜添夜夜添av毛片| 日日撸夜夜添| 久久久久久人妻| 亚洲精品日韩在线中文字幕| 最近手机中文字幕大全| 性高湖久久久久久久久免费观看| 国产亚洲最大av| 丝瓜视频免费看黄片| freevideosex欧美| 亚洲国产色片| 欧美三级亚洲精品| 一区二区三区免费毛片| 女人久久www免费人成看片| 91久久精品国产一区二区三区| av在线app专区| 久久99精品国语久久久| a 毛片基地| 黑丝袜美女国产一区| 成人18禁高潮啪啪吃奶动态图 | 国产精品一区二区在线不卡| 久久久a久久爽久久v久久| 少妇人妻一区二区三区视频| 国产精品99久久99久久久不卡 | 女人久久www免费人成看片| 国产视频首页在线观看| 久久这里有精品视频免费| 一区二区av电影网| 欧美极品一区二区三区四区| av在线观看视频网站免费| 国产亚洲精品久久久com| 女人久久www免费人成看片| 国产高清国产精品国产三级 | 我要看日韩黄色一级片| 联通29元200g的流量卡| 亚洲一级一片aⅴ在线观看| 久久久久久久大尺度免费视频| 91久久精品电影网| av不卡在线播放| 大香蕉久久网| 久久97久久精品| 97超碰精品成人国产| 大码成人一级视频| 又爽又黄a免费视频| 日本午夜av视频| 99热6这里只有精品| 精品熟女少妇av免费看| 日韩免费高清中文字幕av| 黄色配什么色好看| 久久青草综合色| 人妻 亚洲 视频| 黑人高潮一二区| 免费黄网站久久成人精品| 国产精品偷伦视频观看了| 丝袜脚勾引网站| 高清午夜精品一区二区三区| 91精品国产国语对白视频| 国产毛片在线视频| 激情 狠狠 欧美| 国产淫语在线视频| 国产国拍精品亚洲av在线观看| 日日啪夜夜爽| 亚洲综合色惰| 97在线人人人人妻| 少妇人妻久久综合中文| 男人狂女人下面高潮的视频| 久久久午夜欧美精品| 久久久色成人| 三级国产精品欧美在线观看| 午夜老司机福利剧场| 亚洲中文av在线| 青青草视频在线视频观看| 国产在线男女| 国产亚洲午夜精品一区二区久久| 国产免费又黄又爽又色| av免费在线看不卡| 亚洲人成网站在线播| 少妇人妻一区二区三区视频| 久久99热6这里只有精品| 久久热精品热| 欧美精品亚洲一区二区| 国内少妇人妻偷人精品xxx网站| 大香蕉97超碰在线| 免费观看av网站的网址| 美女cb高潮喷水在线观看| h日本视频在线播放| 观看美女的网站| 免费看不卡的av| 亚洲欧洲日产国产| 欧美变态另类bdsm刘玥| 日韩强制内射视频| 欧美日本视频| 18禁在线播放成人免费| 伦精品一区二区三区| 搡老乐熟女国产| 免费人成在线观看视频色| 国产精品人妻久久久影院| 亚洲国产精品国产精品| 日本色播在线视频| 午夜精品国产一区二区电影| 少妇被粗大猛烈的视频| 成年av动漫网址| av免费观看日本| 欧美丝袜亚洲另类| 亚洲精品自拍成人| 一级黄片播放器| 日韩av不卡免费在线播放| av女优亚洲男人天堂| 亚洲av日韩在线播放| 春色校园在线视频观看| 亚洲欧美一区二区三区黑人 | 在线观看一区二区三区激情| 一本—道久久a久久精品蜜桃钙片| 欧美高清性xxxxhd video| 日本午夜av视频| 国产精品伦人一区二区| 日本午夜av视频| 亚洲精品456在线播放app| av播播在线观看一区| a级毛色黄片| 2022亚洲国产成人精品| av卡一久久| 人人妻人人看人人澡| 亚洲欧洲国产日韩| 99热6这里只有精品| 一级毛片黄色毛片免费观看视频| 免费观看性生交大片5| 尾随美女入室| 啦啦啦中文免费视频观看日本| 在线精品无人区一区二区三 | 久久久久久久久久人人人人人人| av国产久精品久网站免费入址| 超碰av人人做人人爽久久| 欧美高清成人免费视频www| av一本久久久久| 国产av一区二区精品久久 | 乱码一卡2卡4卡精品| 精品一区二区免费观看| tube8黄色片| 日韩三级伦理在线观看| 亚洲丝袜综合中文字幕| 亚洲精品乱久久久久久| 亚洲精品中文字幕在线视频 | 免费看光身美女| 国产亚洲午夜精品一区二区久久| 国产在线男女| 国产亚洲最大av| 亚洲欧美成人精品一区二区| 久久人人爽av亚洲精品天堂 | 亚洲欧美清纯卡通| 亚洲欧美精品自产自拍| 亚洲国产精品成人久久小说| 中文字幕av成人在线电影| 最近2019中文字幕mv第一页| 91午夜精品亚洲一区二区三区| 亚洲精品第二区| 免费高清在线观看视频在线观看| 热re99久久精品国产66热6| 少妇人妻一区二区三区视频| 亚洲人成网站在线播| 亚洲国产av新网站| 亚洲av成人精品一区久久| 欧美精品国产亚洲| 成人午夜精彩视频在线观看| 大片免费播放器 马上看| 97热精品久久久久久| 日韩亚洲欧美综合| 欧美日韩精品成人综合77777| 日韩av在线免费看完整版不卡| 黄片无遮挡物在线观看| 国产国拍精品亚洲av在线观看| av一本久久久久| 午夜福利视频精品| 香蕉精品网在线| 99精国产麻豆久久婷婷| 91久久精品国产一区二区三区| 欧美成人精品欧美一级黄| 国产av精品麻豆| 少妇的逼好多水| 香蕉精品网在线| 99精国产麻豆久久婷婷| 99久久精品热视频| 欧美少妇被猛烈插入视频| 春色校园在线视频观看| 韩国高清视频一区二区三区| 人人妻人人澡人人爽人人夜夜| 一边亲一边摸免费视频| 最近的中文字幕免费完整| 色5月婷婷丁香| 在线观看一区二区三区| 熟女人妻精品中文字幕| 这个男人来自地球电影免费观看 | 少妇熟女欧美另类| 国产精品国产三级专区第一集| 国产乱来视频区| 只有这里有精品99| 日本欧美国产在线视频| 亚洲,欧美,日韩| 国产欧美另类精品又又久久亚洲欧美| 秋霞在线观看毛片| h视频一区二区三区| 国产成人freesex在线| 国产黄色免费在线视频| 欧美国产精品一级二级三级 | 日韩欧美精品免费久久| 国产在线男女| 激情五月婷婷亚洲| 看免费成人av毛片| av国产久精品久网站免费入址| 插逼视频在线观看| 日本欧美国产在线视频| 色网站视频免费| 久久精品久久精品一区二区三区| 少妇的逼水好多| 日韩中文字幕视频在线看片 | 国产有黄有色有爽视频| 日本wwww免费看| 最后的刺客免费高清国语| av线在线观看网站| 中文天堂在线官网| 丰满少妇做爰视频| 一级毛片aaaaaa免费看小| 少妇人妻一区二区三区视频| 久久久色成人| .国产精品久久| 国内精品宾馆在线| 亚洲av中文av极速乱| 欧美3d第一页| 99久久精品热视频| 秋霞伦理黄片| av在线播放精品| 色网站视频免费| 国产精品久久久久久精品电影小说 | 在线看a的网站| 亚洲va在线va天堂va国产| 女人久久www免费人成看片| 一级片'在线观看视频| 国产乱来视频区| 国产精品.久久久| 国产真实伦视频高清在线观看| 啦啦啦视频在线资源免费观看| 国产日韩欧美亚洲二区| 国产亚洲5aaaaa淫片| 下体分泌物呈黄色| 女性被躁到高潮视频| 在线观看av片永久免费下载| 啦啦啦啦在线视频资源| 视频区图区小说| 国产色爽女视频免费观看| 麻豆国产97在线/欧美| 尾随美女入室| 美女内射精品一级片tv| 纵有疾风起免费观看全集完整版| 成年免费大片在线观看| 欧美bdsm另类| 国产精品久久久久久久久免| 2022亚洲国产成人精品| 亚洲欧美清纯卡通| 看免费成人av毛片| 麻豆国产97在线/欧美| 中文字幕免费在线视频6| 成人亚洲欧美一区二区av| 多毛熟女@视频| 中文字幕免费在线视频6| 午夜福利在线在线| 日韩 亚洲 欧美在线| 久久久久久九九精品二区国产| 日日摸夜夜添夜夜爱| 亚洲电影在线观看av| 久久久色成人| 亚洲天堂av无毛| 尤物成人国产欧美一区二区三区| 男女边吃奶边做爰视频| 成人高潮视频无遮挡免费网站| 午夜福利影视在线免费观看| 亚洲精品视频女| 亚洲国产成人一精品久久久| 日韩中文字幕视频在线看片 | 免费久久久久久久精品成人欧美视频 | 日本vs欧美在线观看视频 | 女人十人毛片免费观看3o分钟| 国产色爽女视频免费观看| 精华霜和精华液先用哪个| 国产黄色视频一区二区在线观看| 国产精品一区二区在线观看99| 91在线精品国自产拍蜜月| 国产精品人妻久久久久久| 亚洲精品国产色婷婷电影| 五月开心婷婷网| 日韩亚洲欧美综合| 精华霜和精华液先用哪个| 日本欧美视频一区| 国产综合精华液| 熟女电影av网| 三级国产精品片| 在线观看免费日韩欧美大片 | 亚洲伊人久久精品综合| 成年免费大片在线观看| 22中文网久久字幕| 免费观看性生交大片5| 中文字幕人妻熟人妻熟丝袜美| 只有这里有精品99| 久久精品夜色国产| 精品国产一区二区三区久久久樱花 | 亚洲av成人精品一二三区| 国产精品99久久99久久久不卡 | 久久精品国产鲁丝片午夜精品| 国产精品成人在线| 男人爽女人下面视频在线观看| 尾随美女入室| 精品午夜福利在线看| 熟妇人妻不卡中文字幕| 亚洲精品自拍成人| 毛片女人毛片| 人妻制服诱惑在线中文字幕| 超碰av人人做人人爽久久| 亚洲精品一二三| 丰满迷人的少妇在线观看| 亚洲,欧美,日韩| 亚洲av免费高清在线观看| 久久99精品国语久久久| videos熟女内射| 搡老乐熟女国产| 大陆偷拍与自拍| 色视频www国产| av又黄又爽大尺度在线免费看| 国产深夜福利视频在线观看| 2022亚洲国产成人精品| 丝袜脚勾引网站| 你懂的网址亚洲精品在线观看| 国产精品一及| 少妇人妻精品综合一区二区| 国产精品.久久久| 青青草视频在线视频观看| 人人妻人人看人人澡| 大片电影免费在线观看免费| 在线观看一区二区三区激情| 国产成人91sexporn| 久久精品国产鲁丝片午夜精品| www.色视频.com| 新久久久久国产一级毛片| 国产成人91sexporn| 久久av网站| 免费观看的影片在线观看| 欧美日韩在线观看h| 亚洲第一av免费看| 久久青草综合色| 纵有疾风起免费观看全集完整版| 欧美成人午夜免费资源| 国产亚洲最大av| 99视频精品全部免费 在线| 男人狂女人下面高潮的视频| 色婷婷av一区二区三区视频| 成人一区二区视频在线观看| 激情 狠狠 欧美| 国产在线视频一区二区| 嫩草影院入口| 人体艺术视频欧美日本| 免费观看的影片在线观看| 高清黄色对白视频在线免费看 | 天堂8中文在线网| 亚洲,欧美,日韩| 国产成人免费无遮挡视频| 五月开心婷婷网| 亚洲精品国产av成人精品| 国内少妇人妻偷人精品xxx网站| 久久99热6这里只有精品| 亚洲欧美日韩无卡精品| 一级爰片在线观看| 少妇猛男粗大的猛烈进出视频| 亚洲精品视频女| 99热这里只有是精品50| 久久韩国三级中文字幕| 日本与韩国留学比较| 尤物成人国产欧美一区二区三区| 国产精品欧美亚洲77777| 日本黄色日本黄色录像| 好男人视频免费观看在线| 日本-黄色视频高清免费观看| 亚洲四区av|