999精品在线视频,手机成人午夜在线视频,久久不卡国产精品无码,中日无码在线观看,成人av手机在线观看,日韩精品亚洲一区中文字幕,亚洲av无码人妻,四虎国产在线观看 ?

增材制造技術中的流固耦合度影響因素

2021-06-08 18:01:30張闊
工業技術創新 2021年1期
關鍵詞:模型

摘? ?要: 對流固界面的耦合程度為增材制造技術帶來的挑戰進行報道。以流體柱體外圍繞無限遠各向同性固體模型為基準,在角頻率—波數域建立流體滿足的聲波方程、固體滿足的彈性波方程和流固耦合邊界條件。引入流固耦合度這一概念,導出法向位移斷續、法向應力斷續、切向應力為零的流固界面不完全耦合邊界條件,考察不同中心頻率的單極子聲源在模型中激發的頻散特征和波形特征。流體柱體與硬性固體相耦合,降低了各階模式的截止頻率,壓制了各階模式的傳播速度,法向位移耦合度或法向應力耦合度的增大,分別使得各階模式更加遠離或靠近流體柱體外圍不存在硬性固體的情形;流體柱體與軟性固體相耦合,頻散波僅存在一個具有較低截止頻率的0階模式,法向位移耦合度或法向應力耦合度的減小,均使得0階模式更加遠離流體柱體外圍不存在軟性固體的情形,且截止頻率持續降低直至消失,甚至可能導致0階模式發生反轉。法向應力耦合度為零時的復雜波串并非頻散曲線中的某階模式,而是某種其他非頻散波。留存于流體柱體內的“滯能縱波”和近乎于流體柱體外圍繞橫向各向同性固體時的傳播特性,是流體柱體與軟性固體相耦合時產生的特殊數學物理現象。

關鍵詞: 增材制造;流固耦合度;角頻率—波數域;流固界面不完全耦合邊界條件;單極子聲源;頻散特征;截止頻率;0階模式;反轉;非頻散波;“滯能縱波”;橫向各向同性

引言

聲波、彈性波的傳播特性,是評價增材制造工藝效果的基礎指標,是檢驗增材制造技術制備所得物品性能的得力工具[1]。將聲波、彈性波的傳播特性研究貫穿于增材制造全生命周期始終,將有力推動增材制造行業精益生產,促進增材制造產品質量和競爭力提升[2]。

2017年,我國提出行動計劃,明確要求提升增材制造裝備、核心器件及軟件質量[3],其中數學物理方法和算法是支撐增材制造工藝智能規劃、在線檢測、過程控制等軟件功能的核心[4]。聲波、彈性波傳播的數學物理研究方法涵蓋理論法和數值法兩種,理論法包括基于并矢格林函數的求解方法、角頻率—波數域分析法等,數值法包括各時空域的有限差分法、有限元法等。文獻[5]采用融合了完全匹配層全吸收邊界的、適用于斜向各向同性(屬于一種各向異性)介質的時域交錯網格有限差分法,對彈性材料的各向異性為增材制造技術帶來的挑戰進行了報道,呼吁科研人員重視各向異性介質的Thomsen參數大小以及介質在各個方位的傳播速度差異程度對彈性波傳播特性的影響;文獻[6]采用彈性波傳播的并矢格林函數解,以增材制造工藝中常見的負泊松比超材料為研究對象,初步探索了超材料彈性參數線性化的普適方法,為增材制造領域共性技術創新提供了驅動力;文獻[7]采用頻域有限元方法,以1.5維聲波、彈性波方程為例,對內在機理相對復雜、細節把控要求較高、求解具有一定技術門檻的半整數維波動方程進行了考察,尤其是對完全匹配層拉伸尺度選取、等效積分弱形式推導、流固耦合邊界條件設置、極點規避等技術細節進行了注記,為更高維空間波傳播特性的研究打開了局面;文獻[8-9]基于聲波和彈性波在流固界面的雙向轉換機制,推導了應用于時域有限元方法的,帶有完全匹配層的聲波—彈性波耦合(即流固耦合)方程的等效積分弱形式,并立足于增材制造技術應用需求,構思了若干聲波—彈性波耦合模型,研討了遠場、近場條件下聲波—彈性波耦合對波場傳播的影響,得到了大量具有實際價值的結論和啟示。文獻[10]為流固耦合問題提出了一種來源于經驗的傳遞矩陣校正方法,為流固耦合模型頻散曲線的智能求取提供了便利。上述文獻雖然覆蓋了各維度、各時空域下基于各種數學物理方法和算法的聲波、彈性波及二者之間相互耦合機制的研究,但并未考慮到聲波、彈性波之間的耦合程度對波場傳播的影響,而增材制造工藝的固有特點又是導致這種耦合程度的不確定性的根源所在。據此,本文引入流固耦合度這一概念,在角頻率—波數域考察不同流固耦合度下波場傳播的頻散特征和波形特征,以期為增材制造工藝效果評價和物品制備性能檢驗提供更加嚴密、嚴謹的理論基礎。

1? 流體柱體中聲波的傳播

式(40)稱作流固界面不完全耦合邊界條件。其中,、分別稱作法向位移耦合度、法向應力耦合度,處于區間[0, 1];稱作切向應力耦合度,僅取0或1。當、、均為0時,式(40)退化為第1章所討論的流體柱體中聲波的傳播情形;當三者均為1時,式(40)退化為第2章所討論的流體柱體外圍繞無限遠固體,且流固界面完全耦合時,流體柱體內聲波的傳播情形;對于、(二者可以不相等)且的情形,本文將邊界條件描述為:法向位移斷續、法向應力斷續、切向應力為零。

4? 應用與討論

4.1? 模型與參數

將第1章所述流體柱體的情形記作模型1,將第2章所述流體柱體外圍繞無限遠固體的情形記作模型2。限于篇幅,本文僅對單極子聲源激發的波場進行考察。多極子聲源情形的研究結論必然有所不同,但研究思路類似。

模型基本參數如表1所示,其中硬性固體指橫波速度大于相鄰流體中聲波傳播速度的固體,軟性固體反之。

4.2? 頻散特征計算與討論

所有滿足式(46)的點集構成模型2的頻散曲線。對于式(46),可根據文獻[10],首先在給定的波源頻率和軸向相速度范圍內對二維譜進行大體觀察,判斷二維譜的極性翻轉特性,然后針對不同表現形式的極性翻轉特性,采用校正因子實施“反翻轉”,最后根據“反翻轉”過程整理得到傳遞矩陣校正圖譜,從而通過牛頓迭代法[11]求解。牛頓迭代法詳見附錄B。

4.2.1? 流固界面完全耦合

首先,在模型2的流體柱體外圍選用硬性固體,且將流固界面設定為完全耦合(即),繪制頻散曲線,并與模型1相對比,如圖1a所示。對比得知,在0~25 kHz頻段,兩個模型均具有0階、1階、2階模式,其中0階模式無截止頻率,其他階模式具有較高的截止頻率。此外,流體柱體與硬性固體相耦合,不僅降低了各階模式的截止頻率,而且壓制了各階模式的傳播速度。在模型2中,1階以上(含)模式在截止頻率下傳播速度最大,等于硬性固體的橫波速度,而在模型1中,1階以上(含)模式在截止頻率下傳播速度無限大(事實上,模型1的流體柱體模型是一個典型的“剛性壁”波導問題,其等效于流體柱體外圍繞無限遠“波速無限大的超硬性固體”的情形)。

其次,在模型2的流體柱體外圍選用軟性固體,仍將流固界面設定為完全耦合,繪制頻散曲線,并與模型1相對比,如圖1b所示。與圖1a截然不同的是,模型2的頻散曲線似乎與模型1毫無相關性。在0~25 kHz頻段對模型2的頻散曲線進行求解,僅可得到0階模式,且該模式具有一個較低的截止頻率(<1 kHz)。0階模式在截止頻率下傳播速度最大,等于軟性固體的橫波速度。

4.2.2? 流固界面不完全耦合

令人出乎意料的是,對于任意滿足式(47)的流固耦合度,無論模型2采用硬性固體還是軟性固體,計算得到的頻散曲線均與流固界面完全耦合的情形重合。

進一步,考察法向位移耦合度、法向應力耦合度不同的情形。不失一般性,以下僅討論0階、1階模式,并且為便于對比,模型1的上述模式也將在頻散特征分析圖中用虛線繪出。

(1)硬性固體

① 法向位移耦合度不變,法向應力耦合度變化

令不變,將分別取值為0.25、0.5、0.75、1,計算結果如圖2a所示。在的不同取值下,計算結果完全相同,其中模型1和模型2的0階模式發生了重合,1階模式發生了交叉。與圖1a對比得知,頻散曲線交叉現象似乎是由模型2的1階以上(含)模式之間相互交融而引起的。

令不變,將分別取值為0、0.25、0.5、0.75、1,計算結果如圖2b所示。當時,模型2的0階模式不存在,1階模式在截止頻率下傳播速度無限大。隨著的增大,模型2的0階模式開始顯現,并不斷趨近于模型1的0階模式。在的不同取值下,模型2的1階模式近乎巧合地交于一點,當實際頻率一定,且大于這一巧合交點所在的頻率時,傳播速度隨著的增大而增大,反之亦然。初步推斷這一巧合交點所在的頻率為模型1的1階模式的截止頻率,該截止頻率有理論解[12],即

其中,為截止頻率,為第一類變形Bessel函數的第一個虛部大于0的根[13]。將相關參數代入式(48),求得。采用牛頓迭代法在該截止頻率附近對的不同取值進行精算,發現確鑿如此。

分別令、、不變,取值同上,計算結果分別如圖2c、2d、2e所示。除數值必然不同外,上述結果與時的規律相似,其中一個顯著的區別是,當一定時, 的增大,使得模型2的0階模式更加遠離模型1。

② 法向應力耦合度不變,法向位移耦合度變化

限于篇幅,僅討論的情形。將分別取值為0、0.25、0.5、0.75、1,計算結果如圖2f所示。顯而易見,當一定時,的增大,除了使得模型2的0階模式更加遠離模型1以外,也使得模型2的1階模式更加稍微遠離模型1。

(2)軟性固體

章節4.2.1揭示了流體柱體外圍繞無限遠軟性固體時頻散特性的捉摸不定。因此,以下從流固耦合度由大至小開展試探性分析。

① 法向位移耦合度不變,法向應力耦合度變化

令不變,將分別取值為1、0.75、0.5、0.25、0,計算結果如圖3a所示。剛減小到0.75,模型2的0階模式的截止頻率即消失;當減小到0.25時,模型2的0階模式甚至發生了反轉;當時,模型2的0階模式不復存在。

令不變,取值同上,計算結果如圖3b所示。與時不同,當減小到0.5時,模型2的0階模式的截止頻率方消失;和時的頻散特性與時相似。

令不變,取值同上,計算結果如圖3c所示。與以上情形不同,當減小到0.25時,模型2的0階模式的截止頻率方消失,且模型2的0階模式未再發生反轉。

令不變,取值同上,計算結果如圖3d所示。在的不同取值下,模型2的0階模式的截止頻率均未消失,且模型2的0階模式均未發生反轉。此外,模型2的0階模式對的敏感度已變得非常弱。

令不變,將分別取值為1、0.75、0.5、0.25,計算結果如圖3e所示。在的不同取值下,計算結果完全相同,模型2的0階模式頻散極其輕微,且其截止頻率在以上情形中最大(≈3 kHz)。

② 法向應力耦合度不變,法向位移耦合度變化

限于篇幅,僅討論的情形。將分別取值為0.25、0.5、0.75、1,計算結果如圖3f所示。顯而易見,當一定時,的增大,使得模型2的0階模式更加接近流固界面完全耦合的情形。

4.3? 波形特征計算與討論

結合章節4.2計算得到的頻散特征,考察代表性中心頻率3 kHz和10 kHz的單極子聲源在模型1和模型2中激發的波形。聲源與接收器之間的距離取1.0~2.0 m,接收器之間的間隔取0.2 m。

4.3.1? 流固界面完全耦合

(1)硬性固體

繪制中心頻率為3 kHz的單極子聲源在模型1中激發的波形,并與其中直達波的貢獻相對比,如圖4a所示。在模型1中,直達波的貢獻較為微弱,且衰減較快;反射波占據主要地位,為純粹的0階模式。進一步,繪制該聲源在模型2中激發的波形,并與模型1相對比,如圖4b所示。流體柱體外圍硬性固體的存在,使得流體柱體內的模式波發生了遲滯,即傳播速度降低,這與圖1a所示的頻散曲線相符合。

同理,繪制中心頻率為10 kHz的單極子聲源在模型1中激發的波形,并與其中直達波的貢獻相對比,如圖5a所示。在模型1中,直達波的貢獻較為微弱,且衰減更快,但直達波與反射波的幅度比大于聲源中心頻率為3 kHz的情形。盡管如此,反射波依舊占據主要地位,且在一個0階模式脈沖后有復雜波串緊隨,這些波串很容易被認作1階模式。進一步,繪制該聲源在模型2中激發的波形,如圖5b所示。如同章節4.2.1所提到的,流體柱體外圍硬性固體的存在,壓制了各階模式的傳播速度,充其量不會超過硬性固體的橫波速度,這使得模式波特征得到了很大程度的簡化。事實上,圖5b的波串中也蘊含著多種模式,且以非頻散波為主,但這已超出本文的討論范疇,感興趣的讀者可與筆者私下討論。

(2)軟性固體

繪制中心頻率為3 kHz的單極子聲源在模型2中激發的波形,并與模型1相對比,如圖6a所示。流體柱體外圍軟性固體的存在,使得流體柱體內的模式波受到了對沖,幅度嚴重衰減,其他模式波的幅度也遠遠不及流體柱體外圍不存在軟性固體的情形。事實上,模型2的波串中蘊含著兩種模式,一是傳播速度較快的,以固體縱波速度傳播的非頻散波,未在圖1b中展示;二是傳播速度較慢的,以低于固體橫波速度傳播的頻散波,即圖1b中的0階模式。

同理,繪制中心頻率為10 kHz的單極子聲源在模型2中激發的波形,如圖6b所示。此時,波串中仍蘊含著兩種模式,一是幅度較大、傳播速度較快的,以固體縱波速度傳播的非頻散波;二是衰減較嚴重的,與上述非頻散波相交融的,以更加稍微低于固體橫波速度傳播的頻散波。此外,將模型2中的波形與模型1(即圖5a中實線)相對比,發現流體柱體外圍繞軟性固體時產生的非頻散波,似乎是將流體柱體外圍不存在軟性固體時產生的反射波的一部分釋放出去的結果,這使得留存在流體柱體內的這一部分持續了約1 ms左右即趨于平穩。部分文獻將這一非頻散波稱作“漏能縱波”(Leaky P-wave),但按照上述分析,將留存在流體柱體內的這一部分稱作“滯能縱波”(Hysteresis P-wave)似乎更為貼切。

4.3.2? 流固界面不完全耦合

首先考察法向位移耦合度、法向應力耦合度相同的情形。與章節4.2.2的分析相一致,對于任意滿足式(47)的流固耦合度,無論模型2采用硬性固體還是軟性固體,計算得到的波形均與流固界面完全耦合的情形重合。進一步,考察法向位移耦合度、法向應力耦合度不同的情形。

(1)硬性固體

令不變,將分別取值為0.5、1,繪制中心頻率為3 kHz的單極子聲源激發的波形,計算結果與模型1完全相同。進一步,將單極子聲源中心頻率調整至10 kHz,依舊如是。

令不變,將分別取值為0、0.5、1,繪制中心頻率為3 kHz的單極子聲源激發的波形,計算結果對比圖如圖7a和7b所示。當時,復雜波串再次顯現,根據圖2c,3 kHz的中心頻率遠遠小于模型2的1階模式的截止頻率,排除了這一復雜波串為1階模式的可能,它可能是0階模式或某種其他非頻散波。隨著的增大,上述復雜波串消失,此外流體柱體內的模式波傳播速度相應提高,這與圖2c所示的頻散曲線相符合。進一步,將單極子聲源中心頻率調整至10 kHz,計算結果對比圖如圖7c和7d所示。當時,復雜波串再次顯現,但由于模型2的0階模式并不存在,進一步排除了這一復雜波串為0階模式的可能,它只能是某種其他非頻散波。與聲源中心頻率為3 kHz的情形相似,隨著的增大,上述復雜波串消失,此外流體柱體內的模式波傳播速度相應提高,這也與圖2c所示的頻散曲線相符合。但可意外發現,當時,波串持續時間相對更長。

令不變,將分別取值為0、0.5、1,繪制中心頻率為3 kHz、10 kHz的單極子聲源激發的波形,分析結論與時基本相同。限于篇幅,計算結果不再展示。

(2)軟性固體

令不變,將分別取值為0.5、1,繪制中心頻率為3 kHz的單極子聲源激發的波形,計算結果與模型1完全相同。進一步,將單極子聲源中心頻率調整至10 kHz,依舊如是。

令不變,將分別取值為0、0.5、1,繪制中心頻率為3 kHz的單極子聲源激發的波形,計算結果對比圖如圖8a和8b所示。當時,復雜波串再次顯現,只是比硬性固體的情形更規則有序一些,同樣根據圖3c,0階模式并不存在,因此這一復雜波串只能是某種其他非頻散波。隨著的增大,上述復雜波串消失,此外流體柱體內的模式波傳播速度相應提高,這與圖3c所示的頻散曲線相符合,但可以驚奇地發現,、時的模式波竟至于體現為近乎于流體柱體外圍繞無限遠橫向各向同性(屬于一種各向異性)固體時流體柱體內聲波的傳播特性!這極可能造成對流體柱體外圍固體各向異性程度的誤判。進一步,將單極子聲源中心頻率調整至10 kHz,計算結果對比圖如圖8c和8d所示。當時,復雜波串再次顯現,與前文分析相同,其只能是某種其他非頻散波。與聲源中心頻率為3 kHz的情形相似,隨著的增大,上述復雜波串消失,此外流體柱體內的模式波傳播速度相應提高,這也與圖3c所示的頻散曲線相符合,且仍可發現,當時,波串持續時間相對更長。

令不變,將分別取值為0、0.5、1,繪制中心頻率為3 kHz、10 kHz的單極子聲源激發的波形,分析結論與時基本相同。限于篇幅,計算結果不再展示。

5? 結論與局限

本文以流體柱體外圍繞無限遠固體模型為基準,引入流固耦合度這一概念,在角頻率—波數域考察了不同流固耦合度下波場傳播的頻散特征和波形特征,為增材制造工藝效果評價和物品制備性能檢驗提供了更加嚴密、嚴謹的理論基礎。主要結論有:

(1) 流體柱體與硬性固體相耦合,不僅降低了各階模式的截止頻率,而且壓制了各階模式的傳播速度。法向位移耦合度的增大,使得各階模式更加遠離流體柱體外圍不存在硬性固體的情形;法向應力耦合度的增大,使得各階模式更加靠近流體柱體外圍不存在硬性固體的情形。

(2) 流體柱體與軟性固體相耦合,頻散波僅存在0階模式,且該模式具有一個較低的截止頻率。法向位移耦合度或法向應力耦合度的減小,均使得0階模式更加遠離流體柱體外圍不存在軟性固體的情形,且截止頻率持續降低直至消失,甚至可能導致0階模式發生反轉。

(3) 流體柱體無論是與硬性固體還是軟性固體相耦合,法向應力耦合度為0時的復雜波串均并非頻散曲線中的某階模式,而是某種其他非頻散波。

(4) 留存于流體柱體內的“滯能縱波”和近乎于流體柱體外圍繞橫向各向同性固體時的傳播特性,是流體柱體與軟性固體相耦合時產生的特殊數學物理現象。

本研究盡管力求全面,但限于篇幅,仍存在一定的局限,包括但不限于:

(1) 當法向位移耦合度和法向應力耦合度相同且非零時,頻散特征和波形特征并未顯現出差別,表明流固耦合度的完備性有待加強;

(2) 多極子聲源,如偶極子、四極子聲源在流固耦合模型中的傳播特性有待分析;

(3) 流體與各向異性固體相耦合時的波場傳播特性有待考察;

(4) 固體—固體界面的彈性耦合度勢必比流固耦合度更為復雜,有待深入;

(5) 為了有效檢驗由增材制造技術制備所得的更為復雜的物品的性能,流固耦合度、彈性耦合度在數值法,如各時空域的有限差分法、有限元法中的數學物理表達,有待探究。

參考文獻

[1] 邵長金, 楊振清, 周廣剛, 等. 場與波[M]. 東營: 中國石油大學出版社, 2015.

[2] 彼得 S. 潘迪, 羅伯特 P. 紐曼. 六西格瑪管理法: 世界頂級企業追求卓越之道: 第2版[M]. 北京: 機械工業出版社, 2017.

[3] 增材制造產業發展行動計劃(2017-2020年)[R].

[4] 顧樵. 數學物理方法[M]. 北京: 科學出版社, 2019.

[5] 張闊, 劉鶴. 增材制造技術中的彈性各向異性影響因素[J]. 工業技術創新, 2017, 4(4): 53-58.

[6] 張闊. 各向同性負泊松比超材料彈性參數線性化初探[J]. 工業技術創新, 2018, 5(4): 53-60.

[7] 張闊. 關于半整數維波動方程求解的幾點注記[J]. 新一代信息技術, 2019, 2(21): 13-21.

[8] 張闊. 增材制造技術中基于時域有限元方法的聲波—彈性波耦合(一): 理論[J]. 工業技術創新, 2019, 6(5): 74-85, 90.

[9] 張闊. 增材制造技術中基于時域有限元方法的聲波—彈性波耦合(二): 應用[J]. 工業技術創新, 2019, 6(6): 75-86.

[10] 張闊. 柱狀模型傳遞矩陣校正圖譜[J]. 新一代信息技術, 2020, 3(2): 1-9.

[11] 李信真, 車剛明, 歐陽潔, 等. 計算方法: 第2版[M]. 西安: 西北工業大學出版社, 2010.

[12] 張海瀾, 王秀明, 張碧星. 井孔的聲場和波[M]. 北京: 科學出版社, 2004.

[13] 王元明. 數學物理方程與特殊函數[M]. 北京: 高等教育出版社, 2004.

猜你喜歡
模型
一半模型
一種去中心化的域名服務本地化模型
適用于BDS-3 PPP的隨機模型
提煉模型 突破難點
函數模型及應用
p150Glued在帕金森病模型中的表達及分布
函數模型及應用
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
3D打印中的模型分割與打包
主站蜘蛛池模板: 国产精品一老牛影视频| 亚洲精品无码日韩国产不卡| 日韩高清一区 | 国产成人1024精品| 黄片一区二区三区| 国产在线一二三区| 国产h视频免费观看| 日韩在线中文| 国产欧美高清| 午夜视频在线观看区二区| 国产高清毛片| 婷婷开心中文字幕| 强奷白丝美女在线观看| 亚洲欧美日本国产专区一区| 国产偷倩视频| 日韩高清无码免费| 亚洲啪啪网| 人妻91无码色偷偷色噜噜噜| 91国内视频在线观看| 婷婷综合在线观看丁香| 亚洲男人的天堂在线| 91成人精品视频| 波多野结衣无码视频在线观看| 国产无码高清视频不卡| 亚洲欧州色色免费AV| 国产免费黄| 亚洲人成网站在线播放2019| 波多野结衣一二三| 国产精品女主播| 最新亚洲av女人的天堂| 狠狠色成人综合首页| 超清无码一区二区三区| 欧美一级特黄aaaaaa在线看片| 国产欧美高清| 亚洲成人网在线播放| 亚洲精品无码人妻无码| 亚洲一区毛片| 四虎在线观看视频高清无码| 国产久操视频| 制服丝袜 91视频| 又猛又黄又爽无遮挡的视频网站| 在线看AV天堂| 青青草原国产av福利网站| 免费99精品国产自在现线| 久久精品人人做人人综合试看| 亚洲毛片在线看| 最新国产你懂的在线网址| 91小视频在线| 国产在线精彩视频二区| 久久精品人人做人人综合试看| 国产人成在线观看| 无码在线激情片| 在线免费亚洲无码视频| 丁香五月婷婷激情基地| 亚洲人成网站18禁动漫无码| 91日本在线观看亚洲精品| 无码国内精品人妻少妇蜜桃视频| 日韩A∨精品日韩精品无码| 伊在人亚洲香蕉精品播放| 91小视频在线播放| 亚洲国产成人久久77| 欧美.成人.综合在线| 亚洲第一香蕉视频| 国产亚洲视频中文字幕视频| 免费一级无码在线网站| 亚洲av无码成人专区| 亚洲区欧美区| 免费视频在线2021入口| 国产亚洲精品自在久久不卡| 欧美福利在线观看| 香蕉伊思人视频| 综合色天天| 亚洲系列中文字幕一区二区| 亚洲日韩久久综合中文字幕| 亚洲美女AV免费一区| 日韩 欧美 国产 精品 综合| 国产波多野结衣中文在线播放| 成人无码一区二区三区视频在线观看| 在线看片中文字幕| 日本免费一级视频| 免费99精品国产自在现线| 精品剧情v国产在线观看|