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

雙模式熱管堆非核原理樣機(jī)熱力特性數(shù)值分析

2023-08-29 13:00:40張凌義李潘瀟王成龍田文喜蘇光輝秋穗正
原子能科學(xué)技術(shù) 2023年8期

張凌義,李潘瀟,張 胤,王成龍,田文喜,蘇光輝,秋穗正

(西安交通大學(xué) 核科學(xué)與技術(shù)學(xué)院,陜西 西安 710049)

熱管堆的概念在20世紀(jì)60年代由洛斯阿拉莫斯國(guó)家實(shí)驗(yàn)室(LANL)首先提出,其優(yōu)勢(shì)為:系統(tǒng)結(jié)構(gòu)簡(jiǎn)單緊湊;熱管間不存在功能關(guān)聯(lián),安全冗余性好;模塊化程度高,可適配多種熱電熱力轉(zhuǎn)換方式以滿足應(yīng)用場(chǎng)景需求等,具有廣闊的研究前景。美國(guó)Kilopower計(jì)劃提出了千瓦級(jí)小型反應(yīng)堆系統(tǒng)設(shè)計(jì)方案,由高溫?zé)峁茌敓?斯特林裝置發(fā)電,轉(zhuǎn)化效率可達(dá)25%,并于2018年完成了該設(shè)計(jì)的地面演示試驗(yàn)(KRUSTY)[1],目前技術(shù)成熟度已達(dá)到TRL6。國(guó)內(nèi)對(duì)熱管堆的研究相對(duì)處于起步階段,西安交通大學(xué)NuThel研究室于2018年提出了靜默式熱管反應(yīng)堆NUSTER設(shè)計(jì)方案[2-3],并建立了非核原理樣機(jī)[4-5]。李潘瀟等[6]設(shè)計(jì)搭建的雙模式熱管堆非核原理樣機(jī)NUSTER2.0采用“模擬堆芯-高溫?zé)峁?斯特林-溫差發(fā)電”一體化集成設(shè)計(jì),具備高功率與低功率兩種運(yùn)行模式,可適配多樣化能源需求。本文基于COMSOL以NUSTER2.0原理樣機(jī)為對(duì)象搭建系統(tǒng)耦合仿真模型,開(kāi)展穩(wěn)態(tài)實(shí)驗(yàn)驗(yàn)證,并對(duì)其運(yùn)行特性及力學(xué)特性進(jìn)行仿真研究。

1 原理樣機(jī)簡(jiǎn)介及穩(wěn)態(tài)實(shí)驗(yàn)研究

雙模式熱管堆非核原理樣機(jī)NUSTER2.0由固態(tài)堆芯模塊、高溫?zé)峁苣K、動(dòng)態(tài)熱聲發(fā)電模塊、靜態(tài)溫差發(fā)電模塊、保護(hù)氣腔、余熱排出模塊、控制模塊及數(shù)據(jù)采集模塊組成,樣機(jī)概念圖及布置示意圖如圖1、2所示。該設(shè)計(jì)采用了動(dòng)靜結(jié)合的熱電轉(zhuǎn)換方式,低功率模式下,由溫差發(fā)電靜態(tài)模塊執(zhí)行熱電轉(zhuǎn)換,運(yùn)行噪聲小;高功率模式下,靜態(tài)模塊與斯特林動(dòng)態(tài)模塊共同進(jìn)行發(fā)電,轉(zhuǎn)換效率高。

圖1 非核原理樣機(jī)概念圖

圖2 非核原理樣機(jī)布置示意圖

加熱棒及高溫?zé)峁馨惭b在堆芯基體的孔隙中,加熱棒與堆芯基體間涂有導(dǎo)熱膠以強(qiáng)化導(dǎo)熱。熱量通過(guò)37根高溫鉀金屬?gòu)澱蹮峁茌斔偷届o態(tài)模塊集熱基體及斯特林熱端進(jìn)行熱電轉(zhuǎn)換。靜態(tài)模塊包括集熱基體、5×6布置的銻化鉍溫差發(fā)電片及調(diào)節(jié)氣板,動(dòng)態(tài)模塊為自由活塞式斯特林裝置。基于LabVIEW搭建的測(cè)試平臺(tái)可以對(duì)熱電偶、電子負(fù)載儀和功率儀測(cè)得的數(shù)據(jù)進(jìn)行實(shí)時(shí)的采集與記錄。樣機(jī)15 kW加熱功率穩(wěn)態(tài)實(shí)驗(yàn)結(jié)果顯示,原理樣機(jī)基體平均溫度為627 ℃,最高溫度為635 ℃,低于系統(tǒng)設(shè)計(jì)安全限值850 ℃;氣板平均溫度為260.2 ℃,斯特林基體溫度為475.3 ℃;動(dòng)、靜態(tài)模塊熱管蒸發(fā)段到冷凝段溫差分別為74 K與120 K;動(dòng)態(tài)模塊發(fā)電功率為2 019 W,靜態(tài)模塊發(fā)電功率為493.2 W,整體熱電轉(zhuǎn)換效率為16.7%;動(dòng)、靜態(tài)模塊運(yùn)行參數(shù)實(shí)驗(yàn)結(jié)果列于表1、2。本實(shí)驗(yàn)驗(yàn)證了樣機(jī)設(shè)計(jì)的合理性,并為仿真模型提供了比對(duì)數(shù)據(jù)。

表1 靜態(tài)模塊15 kW穩(wěn)態(tài)實(shí)驗(yàn)結(jié)果

表2 動(dòng)態(tài)模塊15 kW穩(wěn)態(tài)實(shí)驗(yàn)結(jié)果

2 系統(tǒng)耦合仿真模型構(gòu)建及驗(yàn)證

2.1 數(shù)學(xué)物理模型

1) 高溫?zé)峁軣嶙杈W(wǎng)絡(luò)模型

高溫?zé)峁芡ǔ7譃楣鼙趨^(qū)、吸液芯、蒸汽區(qū)3個(gè)區(qū)域,本文通過(guò)如圖3所示的熱阻網(wǎng)絡(luò)模型對(duì)其總體傳熱過(guò)程進(jìn)行描述。可看出,管壁僅受導(dǎo)熱方程控制,而吸液芯與蒸氣區(qū)傳熱機(jī)理較為復(fù)雜。

圖3 熱管熱阻網(wǎng)絡(luò)示意圖

在吸液芯區(qū),假定其工質(zhì)完全熔化,忽略其內(nèi)部液體的流動(dòng),將工質(zhì)和吸液芯整體簡(jiǎn)化為純導(dǎo)熱模型[7],等效參數(shù)由Chi模型[8-9]求解:

(1)

(2)

其中:ρeff為吸液芯等效密度;ceff為吸液芯等效比熱容;keff為吸液芯等效導(dǎo)熱系數(shù);ε為吸液芯孔隙率;ρ1為液態(tài)鉀密度;c1為液態(tài)鉀比熱容;k1為液態(tài)鉀導(dǎo)熱系數(shù);ρw為吸液芯固體材料密度;cw為吸液芯固體材料比熱容;kw為吸液芯固體材料導(dǎo)熱系數(shù)。

在蒸汽區(qū),當(dāng)熱管完全啟動(dòng)時(shí)蒸汽為連續(xù)流動(dòng)態(tài),主要包括蒸發(fā)段氣液界面相變傳熱熱阻R1、蒸汽軸向流動(dòng)傳熱熱阻R2和冷凝段氣液界面相變傳熱熱阻R3,其中R2占主導(dǎo):

(3)

其中:TV為蒸汽溫度;r為汽化潛熱;Le為熱管的有效長(zhǎng)度,即蒸發(fā)段和冷凝段的一半長(zhǎng)度與絕熱段的和;μV為蒸汽動(dòng)力黏度系數(shù);ρV為蒸汽密度。

2) 溫差發(fā)電裝置模型

溫差發(fā)電裝置一方面存在熱電效應(yīng),一方面因結(jié)構(gòu)特點(diǎn)存在熱對(duì)流與熱輻射現(xiàn)象,討論其傳熱特性時(shí)需要假設(shè)忽略湯姆遜效應(yīng)、熱對(duì)流現(xiàn)象及熱輻射現(xiàn)象,主要考慮其等效導(dǎo)熱能力,理論模型如圖4所示。

圖4 溫差發(fā)電裝置模型

溫差發(fā)電器包括n個(gè)熱電偶對(duì),因此每個(gè)溫差發(fā)電單元的導(dǎo)熱系數(shù)由n個(gè)熱電偶對(duì)的熱電偶對(duì)等效導(dǎo)熱系數(shù)K1、導(dǎo)流片等效導(dǎo)熱系數(shù)K2、陶瓷基底等效導(dǎo)熱系數(shù)K3及空氣間隙等效導(dǎo)熱系數(shù)K4決定,其計(jì)算方法為:

(4)

(5)

(6)

(7)

其中:kp、kn為單熱電偶對(duì)中熱電偶P極、N極導(dǎo)熱系數(shù);kd為單熱電偶對(duì)中導(dǎo)流銅片導(dǎo)熱系數(shù);k1為單熱電偶對(duì)中陶瓷基板導(dǎo)熱系數(shù);k2為單熱電偶對(duì)中空氣間隙導(dǎo)熱系數(shù);A、A1、A2為單熱電偶對(duì)中熱電偶、導(dǎo)流銅片及陶瓷基板橫截面積;L、L1、L2為單熱電偶對(duì)中熱電偶、導(dǎo)流銅片及陶瓷基板高度。

考慮熱電偶與導(dǎo)流片、導(dǎo)流片與基體、基板與熱端、基板與冷端的接觸熱阻R1、R2、RH、RC,由式(8)可求得熱電片等效導(dǎo)熱系數(shù)。

K=K1+K2+K3+K4+

(8)

3) 斯特林裝置模型

自由活塞式斯特林發(fā)動(dòng)機(jī)內(nèi)部流動(dòng)為可壓縮工質(zhì)的非定常交變流動(dòng),對(duì)其直接模擬計(jì)算非常困難,因此需要對(duì)其內(nèi)部流動(dòng)進(jìn)行簡(jiǎn)化假設(shè),通常采用絕熱分析法[10]。本文對(duì)斯特林裝置進(jìn)行如圖5所示的簡(jiǎn)化建模,下腔室即膨脹腔,上腔室即壓縮腔,并假定其體平均溫度為工作腔溫度。

圖5 斯特林裝置模型

為了描述斯特林在系統(tǒng)傳熱場(chǎng)中的影響,本文在斯特林絕熱分析程序中以100 K為步長(zhǎng),保持工質(zhì)物性等輸入?yún)?shù)恒定,在系統(tǒng)工作溫度范圍內(nèi)遞進(jìn)改變膨脹腔溫度Tc和壓縮腔溫度Te,求得軸功率W關(guān)于Tc、Te的二元擬合函數(shù)式,代入模型整體熱耗率P0中,實(shí)現(xiàn)與斯特林程序間的數(shù)據(jù)耦合。

P0=W(Tc,Te)

300 K

(9)

4) 其他組件及系統(tǒng)邊界條件

除高溫?zé)峁堋犭娖八固亓盅b置外的其他組件,均視為僅受導(dǎo)熱方程控制。系統(tǒng)邊界條件包括組件耦合處的接觸邊界條件、加熱通孔施加的熱通量熱源、熱電轉(zhuǎn)換產(chǎn)生的熱耗以及冷卻模塊的熱阱,其余部分視為絕熱邊界條件。

2.2 仿真結(jié)果驗(yàn)證

將熱源邊界條件設(shè)定為15 kW進(jìn)行穩(wěn)態(tài)仿真研究,在該工況下的系統(tǒng)溫度場(chǎng)如圖6所示,樣機(jī)系統(tǒng)溫度場(chǎng)峰值溫度為890 K,發(fā)生在堆芯基體中心處。仿真與實(shí)驗(yàn)得到的高溫?zé)峁芊植记€如圖7、8所示,系統(tǒng)傳熱角度吻合良好,動(dòng)態(tài)模塊中冷熱端總溫降誤差為14.5 K,靜態(tài)模塊中總溫降誤差為17.1 K。

圖6 15 kW工況溫度場(chǎng)分布

a——仿真;b——實(shí)驗(yàn)

a——仿真;b——實(shí)驗(yàn)

對(duì)關(guān)鍵運(yùn)行參數(shù)進(jìn)行進(jìn)一步比對(duì)分析(表3、4),主要參數(shù)中冷凝段平均溫度與調(diào)節(jié)氣板平均溫度相對(duì)誤差較大,最大達(dá)到6.73%。分析認(rèn)為,前者誤差產(chǎn)生的原因是計(jì)算熱管氣腔導(dǎo)熱系數(shù)時(shí)所取等效長(zhǎng)度等為軸向高度,而實(shí)際彎折熱管軸線長(zhǎng)度大于軸向高度,導(dǎo)致模型的導(dǎo)熱系數(shù)偏小。此外,動(dòng)態(tài)模塊熱管溫度曲線吻合程度相對(duì)于靜態(tài)模塊較低,這是由于斯特林程序與系統(tǒng)模型耦合程度不足,存在一定誤差。后者誤差產(chǎn)生的原因是氣板作為系統(tǒng)中溫度梯度最大的區(qū)域,實(shí)驗(yàn)中由兩面熱電偶測(cè)點(diǎn)得到的溫度與仿真結(jié)果中取的體平均溫度存在儀器誤差。除去系統(tǒng)誤差,主要參數(shù)最大相對(duì)誤差不超過(guò)6.29%,驗(yàn)證了仿真模型的可靠性。

表3 15 kW工況靜態(tài)模塊參數(shù)仿真結(jié)果比對(duì)

表4 15 kW工況動(dòng)態(tài)模塊參數(shù)仿真結(jié)果比對(duì)

3 原理樣機(jī)運(yùn)行特性及力學(xué)特性仿真研究

3.1 變功率穩(wěn)態(tài)仿真研究

原理樣機(jī)各組件的運(yùn)行溫度存在較大區(qū)別,設(shè)計(jì)上通過(guò)集熱基體及調(diào)節(jié)氣板初步實(shí)現(xiàn)適配,但主要是基于局部測(cè)試的結(jié)果進(jìn)行設(shè)計(jì),缺乏系統(tǒng)整體數(shù)值支撐。對(duì)加熱功率使用敏感性分析得到了不同穩(wěn)態(tài)下的仿真結(jié)果,關(guān)鍵運(yùn)行參數(shù)隨加熱功率變化曲線如圖9所示。由于系統(tǒng)以導(dǎo)熱為主,曲線呈現(xiàn)較好的線性關(guān)系。結(jié)合組件運(yùn)行特性分析,最早對(duì)運(yùn)行形成限制的為動(dòng)態(tài)模塊,在16.5 kW加熱功率下,斯特林基體達(dá)到撞桿現(xiàn)象臨界溫度520 ℃即793.15 K,將導(dǎo)致效率降低、部件損壞等后果,須調(diào)整動(dòng)態(tài)模塊水冷回路中冷卻水流量,保證斯特林裝置在安全范圍內(nèi)運(yùn)行。將不同功率下的熱電片冷、熱端溫度Tc與Th代入效率測(cè)試曲線,如圖10所示,可知25 kW附近靜態(tài)模塊熱電效率達(dá)到峰值,可在保證其他組件冷卻充分的前提下保持加熱功率處于該值附近以取得較高熱電效率。

圖9 樣機(jī)組件溫度-功率曲線

圖10 靜態(tài)模塊效率曲線

3.2 系統(tǒng)熱力耦合穩(wěn)態(tài)仿真研究

熱管冷卻反應(yīng)堆因其固態(tài)屬性,在高溫運(yùn)行下由熱膨脹引起的反應(yīng)性負(fù)反饋不可忽略,因此,本文在考慮熱彈性阻尼造成的能量損失前提下,通過(guò)熱膨脹多物理場(chǎng)接口進(jìn)行熱力強(qiáng)耦合研究,探究15 kW加熱功率下堆芯及其他關(guān)鍵部件的應(yīng)力場(chǎng)及位移場(chǎng)。模擬中將加熱棒視為剛性材料,其所在孔隙受固定約束。

堆芯基體及熱電片的熱膨脹形變?nèi)鐖D11所示,堆芯所受最大等效應(yīng)力為160 MPa,發(fā)生在熱管耦合處,堆芯體膨脹率為1.21%,結(jié)合燃料分布可在蒙特卡羅程序中計(jì)算熱膨脹效應(yīng)反饋系數(shù)aE。熱電片所受最大等效應(yīng)力為101 MPa,危險(xiǎn)區(qū)主要在熱端與調(diào)節(jié)氣板結(jié)合處,熱膨脹形變將一定程度影響其內(nèi)部PN結(jié)的連接,其對(duì)熱電片性能的具體影響需要進(jìn)行單獨(dú)研究。

斯特林裝置的von-Mises應(yīng)力分布如圖12所示,危險(xiǎn)區(qū)主要集中在隔熱環(huán)與基體接觸帶上,等效應(yīng)力達(dá)到366 MPa,發(fā)生的原因是隔熱環(huán)為陶瓷材料,尺寸穩(wěn)定性強(qiáng),在實(shí)裝中隔熱環(huán)設(shè)有自調(diào)裕量,將極大緩解該應(yīng)力。熱管與熱電裝置的耦合處也存在較大接觸應(yīng)力,最大為97 MPa,該應(yīng)力實(shí)際上利于耦合固定。

a——斯特林基體;b——熱管-集熱基體耦合處

除隔熱環(huán)外,15 kW加熱功率下各組件等效應(yīng)力峰值均符合von-Mises屈服條件,且具有150 MPa以上安全裕度,滿足強(qiáng)度要求。提高實(shí)驗(yàn)功率前應(yīng)針對(duì)危險(xiǎn)區(qū)提前進(jìn)行仿真評(píng)估,避免組件發(fā)生屈服失效乃至斷裂,引發(fā)嚴(yán)重后果。

3.3 啟動(dòng)瞬態(tài)仿真研究

由于熱管堆的傳熱特性,原理樣機(jī)的啟動(dòng)方案將對(duì)各組件是否實(shí)現(xiàn)平穩(wěn)啟動(dòng)具有較大影響。為探索較合理的啟動(dòng)策略,對(duì)系統(tǒng)整體進(jìn)行了“3-5-10-15 kW”階梯式功率平臺(tái)啟動(dòng)方案的瞬態(tài)仿真研究,該方案即在不同時(shí)間節(jié)點(diǎn)階躍式地提升功率大小以逐步啟動(dòng)樣機(jī)。樣機(jī)啟動(dòng)特性曲線如圖13所示,各熱電轉(zhuǎn)換功能組件在每個(gè)功率階梯都已基本實(shí)現(xiàn)熱響應(yīng),并在其平均溫度幾乎趨于穩(wěn)態(tài)時(shí)進(jìn)入下一階梯,整體上形成了高效率的平滑啟動(dòng)曲線,約16.7 h系統(tǒng)完全啟動(dòng)。

圖13 啟動(dòng)特性曲線

結(jié)合工程實(shí)際,在“10-15 kW”功率階梯之間設(shè)計(jì)了短暫間隙以減輕加熱棒負(fù)荷。總體而言,功率平臺(tái)方案有利于高溫?zé)峁軆?nèi)部工質(zhì)的熔融,可避免加熱過(guò)快導(dǎo)致局部蒸干等失效工況發(fā)生;同時(shí),該方案中斯特林基體溫升速率低于0.01 K/s,根據(jù)前期斯特林獨(dú)立測(cè)試實(shí)驗(yàn)經(jīng)驗(yàn),此溫度變化速率下斯特林內(nèi)部精細(xì)連接件在啟動(dòng)中受熱沖擊松脫的可能性極低,可以有效延長(zhǎng)其穩(wěn)態(tài)使用壽命。在連續(xù)完成多種工況實(shí)驗(yàn)時(shí),保持低功率模式進(jìn)行過(guò)渡,維持熱管啟動(dòng)狀態(tài),可實(shí)現(xiàn)系統(tǒng)整體的快速熱啟動(dòng)。因此對(duì)于雙模式非核原理樣機(jī),功率平臺(tái)冷啟動(dòng)與保溫?zé)釂?dòng)結(jié)合是較理想的啟動(dòng)策略。

4 結(jié)論

本文在COMSOL中建立了雙模式熱管堆原理樣機(jī)NUSTER2.0的系統(tǒng)耦合仿真模型,并開(kāi)展了樣機(jī)15 kW穩(wěn)態(tài)實(shí)驗(yàn)對(duì)仿真結(jié)果進(jìn)行初步驗(yàn)證。針對(duì)具體工況,對(duì)樣機(jī)的運(yùn)行特性和熱力耦合特性進(jìn)行了仿真研究,得到如下結(jié)論。

1) 15 kW穩(wěn)態(tài)工況下仿真結(jié)果與實(shí)驗(yàn)結(jié)果吻合較好,主要運(yùn)行參數(shù)最大相對(duì)誤差為6.73%,在忽略實(shí)驗(yàn)系統(tǒng)誤差的情況下,相對(duì)誤差不超過(guò)6.29%,初步驗(yàn)證了該模型可對(duì)耦合系統(tǒng)整體特性進(jìn)行準(zhǔn)確仿真。

2) 變功率穩(wěn)態(tài)工況確定了樣機(jī)在提升加熱功率過(guò)程中達(dá)到16.5 kW時(shí)斯特林裝置面臨撞桿臨界,在更高功率水平運(yùn)行時(shí),需要通過(guò)提高動(dòng)態(tài)模塊冷卻水流量強(qiáng)化其熱阱;并確定了靜態(tài)模塊在25 kW附近達(dá)到其效率峰值。

3) 啟動(dòng)瞬態(tài)工況得到了系統(tǒng)在“3-5-10-15 kW”功率平臺(tái)下的啟動(dòng)特性曲線,約16.7 h樣機(jī)實(shí)現(xiàn)完全啟動(dòng)。該方案有效保護(hù)了斯特林精細(xì)部件,并使熱管內(nèi)工質(zhì)充分熔融,安全啟動(dòng);歸納得到功率平臺(tái)冷啟動(dòng)與保溫?zé)釂?dòng)結(jié)合是樣機(jī)較理想的啟動(dòng)策略。

4) 熱力耦合穩(wěn)態(tài)工況求得15 kW加熱功率下熱管與堆芯基體耦合處最大等效應(yīng)力為160 MPa,熱膨脹率為1.21%;熱管與集熱基體耦合處最大接觸應(yīng)力為97 MPa,熱電片所受最大等效應(yīng)力為101 MPa,斯特林隔熱環(huán)的自調(diào)裕量緩解了其與斯特林基體間366 MPa的等效應(yīng)力。

主站蜘蛛池模板: 无码内射在线| 中文无码伦av中文字幕| 欧美精品黑人粗大| 亚洲精品无码av中文字幕| 激情无码视频在线看| 国产精品网曝门免费视频| 久久国产精品夜色| m男亚洲一区中文字幕| 国产欧美日韩精品综合在线| 99激情网| 在线观看国产小视频| 一级毛片在线免费视频| 色香蕉网站| 亚洲日韩精品伊甸| 国产伦精品一区二区三区视频优播| 波多野结衣中文字幕一区二区| 亚洲国产精品不卡在线| 日韩美一区二区| 97人妻精品专区久久久久| 久久99精品久久久大学生| 国产成人麻豆精品| 亚洲无码视频喷水| 国产女同自拍视频| 一本大道无码高清| 精久久久久无码区中文字幕| 亚洲第一香蕉视频| 日本成人精品视频| 日韩免费成人| 免费全部高H视频无码无遮掩| 91蝌蚪视频在线观看| 精品亚洲国产成人AV| 国产成a人片在线播放| 久久精品只有这里有| 国产成人一区免费观看| 国产无码精品在线播放| 日本三区视频| 欧美国产日韩在线播放| 露脸一二三区国语对白| 欧美午夜在线视频| 中文精品久久久久国产网址| 国产女人在线观看| 激情综合五月网| 中文字幕久久亚洲一区| a欧美在线| 在线国产91| 亚洲高清在线播放| 国产91久久久久久| 亚洲欧美日韩综合二区三区| 狠狠色丁香婷婷| 久久精品国产91久久综合麻豆自制| 国产浮力第一页永久地址| 亚洲欧美日韩动漫| 国产精品3p视频| 色综合国产| 国产在线精彩视频论坛| 毛片大全免费观看| 中文字幕首页系列人妻| 欧美第一页在线| 国产国语一级毛片| 婷婷色一区二区三区| 国产一二三区视频| 国产jizz| 亚洲国产成人久久77| 色偷偷一区二区三区| A级毛片无码久久精品免费| 亚洲一区网站| 美女被操黄色视频网站| 91精品国产一区| 国产精品极品美女自在线看免费一区二区| 18黑白丝水手服自慰喷水网站| 国产电话自拍伊人| 欧洲亚洲欧美国产日本高清| www.99精品视频在线播放| 国产电话自拍伊人| 中文字幕永久在线看| 日本免费新一区视频| 小说区 亚洲 自拍 另类| 欧美精品导航| 亚洲中文字幕23页在线| 亚洲综合婷婷激情| 免费国产好深啊好涨好硬视频| 巨熟乳波霸若妻中文观看免费 |