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

    高水頭混流式水輪機減負荷瞬態(tài)流動特性研究

    2022-12-05 06:29:00孫龍剛徐卓飛郭鵬程鄭小波
    水利學(xué)報 2022年11期
    關(guān)鍵詞:導(dǎo)葉轉(zhuǎn)輪水輪機

    孫龍剛,徐卓飛,郭鵬程,2,鄭小波,2

    (1.西安理工大學(xué) 水利水電學(xué)院,陜西 西安 710048;2.西安理工大學(xué) 省部共建西北旱區(qū)生態(tài)水利國家重點實驗室,陜西 西安 710048)

    1 研究背景

    水電由于具有快速高效的啟停響應(yīng),因而被用來調(diào)節(jié)電網(wǎng)容量及頻率的動態(tài)平衡。近年來清潔能源在電網(wǎng)中占比持續(xù)提高,以風(fēng)能和太陽能為代表的清潔能源的不定時并網(wǎng)對電網(wǎng)造成劇烈的不穩(wěn)定沖擊,此時水電機組必然經(jīng)歷快速且頻繁的工況轉(zhuǎn)換如減負荷過程以平衡電網(wǎng)參數(shù)[1-2]。水輪機運行模式的轉(zhuǎn)變,通常伴隨著回流、流動分離、空化、尾水管渦帶、葉道渦、卡門渦等不穩(wěn)定流動現(xiàn)象[3-6],從而在機組內(nèi)部誘發(fā)具有繁雜幅頻特征的壓力脈動。高水頭水輪機由于所受動態(tài)液體力更大,且水輪機減負荷過程中轉(zhuǎn)輪對來流的不適應(yīng)性增強,機組將承受更加復(fù)雜的不平衡作用力,由此產(chǎn)生的交變應(yīng)力載荷會進一步加劇機組關(guān)鍵部件的疲勞破壞和磨損,縮短機組壽命,降低機組水力性能,從而影響電站的安全穩(wěn)定運行[7-9]。因此,開展減負荷過程中高水頭水輪機動態(tài)響應(yīng)研究,對進一步促進水電與其他清潔能源的多能互補具有重要的現(xiàn)實意義。

    模型試驗與數(shù)值模擬是明確工況轉(zhuǎn)換過程中水輪機內(nèi)部不穩(wěn)定流動特性及其誘發(fā)水力振動的重要手段。Trivedi等[10]試驗研究了變轉(zhuǎn)速工況下功率急速上升和下降過程中高水頭混流式水輪機內(nèi)部的壓力脈動響應(yīng),發(fā)現(xiàn)動靜干涉頻率對應(yīng)的壓力脈動幅值在負荷變化期間相對較小,而在瞬態(tài)過程結(jié)束時提高了30倍。此外,快速啟停策略下水輪機轉(zhuǎn)輪葉片最大壓力脈動幅值是緩慢啟停策略幅值的1.5倍[11]。周勤等[12]通過數(shù)值方法研究水泵水輪機,發(fā)現(xiàn)甩負荷過程中轉(zhuǎn)輪進口回流是不穩(wěn)定壓力脈動及不均勻湍動能生成的重要原因,葉片扭矩和徑向力的最大幅值分別達初始值的10倍和60倍。Fu等[13]基于三維大渦模擬技術(shù)捕捉到水泵水輪機內(nèi)部的異常低頻壓力脈動,研究發(fā)現(xiàn)該種壓力脈動的誘發(fā)與水輪機內(nèi)部水錘現(xiàn)象、轉(zhuǎn)輪進口局部回流以及無葉區(qū)擋水環(huán)有關(guān)。基于一維管道與三維數(shù)值模擬耦合方法,陳秋華等[14]發(fā)現(xiàn)以部分負荷工況為初始工況,水泵水輪機飛逸過程中產(chǎn)生的轉(zhuǎn)輪徑向力和壓力脈動相對額定工況更為劇烈。Goyal等[15-16]通過水輪模型試驗發(fā)現(xiàn),甩負荷過程中水輪機尾水管內(nèi)形成的螺旋渦帶可分解為表征軸向運動的突進模式和表征軸向運動的旋轉(zhuǎn)模式,且突進模式比旋轉(zhuǎn)模式提前出現(xiàn),而增負荷過程中柱狀渦帶的形成主要受負荷變化過程中尾水管內(nèi)駐點、反向流動及回流區(qū)的直接影響。Hosseinimanesh等[17]的研究顯示,混流式水輪機運行在空載條件下,尾水管內(nèi)部的軸向反向流動和周向切向流動之間交界面上出現(xiàn)的不穩(wěn)定強剪切層,是造成渦結(jié)構(gòu)破裂和壓力波動向上游傳播至轉(zhuǎn)輪葉片的主要原因。Liu等[18]數(shù)值研究結(jié)果表明,水泵水輪機甩負荷過程中,尾水管壁面附近出現(xiàn)大量大尺度渦結(jié)構(gòu),引發(fā)低頻壓力脈動。Liu等[19]對原型水泵水輪機甩負荷條件下渦帶的形成過程進行了分析,發(fā)現(xiàn)尾水管中的渦帶出現(xiàn)在零力矩運行之前,隨著流量的減小,尾水管中的渦帶強度增加。作者前期研究顯示[20-21],水輪機甩負荷過程中,尾水管內(nèi)的渦旋強度顯著提升,而飛逸過程中的能量耗散主要發(fā)生在轉(zhuǎn)輪和尾水管內(nèi)部。此外,F(xiàn)u等[22-24]基于數(shù)值模擬方法研究了間隙、空化以及水聲學(xué)對水泵水輪機工況轉(zhuǎn)換過程的影響。

    綜上所述,水輪機工況轉(zhuǎn)換會顯著提升其內(nèi)部的壓力脈動幅值,且直接影響水輪機尾水管內(nèi)的渦流結(jié)構(gòu)形態(tài)及其演化過程。為進一步明確工況轉(zhuǎn)換條件下水輪機特征參數(shù)的動態(tài)響應(yīng)、水輪機內(nèi)部壓力脈動及渦流演化規(guī)律,本文基于動網(wǎng)格技術(shù)對一高水頭混流式模型水輪機負荷減小過程進行了瞬態(tài)數(shù)值研究,分析了水輪機壓力和速度演變規(guī)律、尾水管內(nèi)渦流演化特性,并建立轉(zhuǎn)輪軸向力與不穩(wěn)定壓力脈動之間的關(guān)聯(lián)。

    2 數(shù)值模擬策略

    2.1 高水頭模型水輪機本文以挪威科技大學(xué)(Norwegian University of Science and Technology-NTNU)公開的Francis99模型水輪機為研究對象,開展水輪機負荷減小過程的瞬態(tài)數(shù)值研究[25-27]。Francis99模型水輪機三維幾何結(jié)構(gòu)如圖1所示,由進口至出口分別為帶有14個嵌入式固定導(dǎo)葉的蝸殼、葉片數(shù)為28的活動導(dǎo)葉、由15個分流葉片和15個長葉片組成的轉(zhuǎn)輪以及彎肘型尾水管。該水輪機模型與原型轉(zhuǎn)輪出口直徑分別為0.349 m和1.778 m,對應(yīng)的原型機和模型機水頭分別為12.0 m和377 m。圖2所示為Francis99模型水輪機負荷減小過程中活動導(dǎo)葉開度隨時間變化曲線,由圖可知,負荷減小過程涉及到3種運行工況,分別為時段t=0至t=1.0 s,活動導(dǎo)葉開度為的9.84°的最優(yōu)工況(Best Efficiency Point-BEP),時段t=3.5至t=10.0 s,導(dǎo)葉開度為6.72°的部分負荷工況(Part Load-PL)以及時段t=1.0至t=3.5 s,導(dǎo)葉開度由9.84°線性減小至6.72°的負荷減小過程(BEP to PL)。此外,整個負荷變化過程中,水輪機轉(zhuǎn)輪的轉(zhuǎn)速均為333.0 r/min。

    圖1 Francis99模型水輪機

    圖2 活動導(dǎo)葉開度變化規(guī)律

    為了對數(shù)值方法及模擬控制進行驗證,NTNU提供了部分模型試驗結(jié)果,其中壓力及速度測量位置如圖3所示[26-27]。圖3中共布置了6個壓力測點,其中VL02位于活動導(dǎo)葉與轉(zhuǎn)輪之間的無葉區(qū),用以捕捉動靜干涉效應(yīng);轉(zhuǎn)輪葉片壓力面中間位置布置測點P42,出口側(cè)靠近上冠處布置測點P71,用來分析轉(zhuǎn)輪內(nèi)部的壓力脈動特征;同理,轉(zhuǎn)輪葉片吸力面靠近上冠處布置測點S51;尾水管錐管段間隔180°布置了兩個壓力測點DT5和DT6,用來分析尾水管壓力場隨負荷變化規(guī)律。此外,試驗中通過PIV(Particle Image Velocimetry)技術(shù)對尾水管錐管段3條線段L1、L2和L3進行了速度測量,其中L1與L2為垂直距離相差120 mm且與水輪機軸線垂直的水平線,其與x軸之間的角度均為79°;L3為與水輪機軸線重合的豎直線。速度測量中,對L1和L2上的28個位置進行測量,而L3上的測點數(shù)為19。

    圖3 壓力與速度測點布置示意圖

    2.2 導(dǎo)葉網(wǎng)格變形控制水輪機工況轉(zhuǎn)換是通過改變活動導(dǎo)葉開度進而控制水輪機流量來實現(xiàn)的,因此水輪機負荷改變數(shù)值研究中的一項重要內(nèi)容為按照一定的導(dǎo)葉運動規(guī)律動態(tài)調(diào)整活動導(dǎo)葉的開度且保證一定的網(wǎng)格精度。為了實現(xiàn)上述目標,本文采用一種幾何結(jié)構(gòu)布爾運算與網(wǎng)格重構(gòu)相結(jié)合的方法為活動導(dǎo)葉模擬域提供網(wǎng)格的自動更新,如圖4和圖5所示。當(dāng)導(dǎo)葉開度變化時,導(dǎo)葉區(qū)域的網(wǎng)格發(fā)生變形,求解器按照預(yù)設(shè)的網(wǎng)格質(zhì)量判斷準則,確定是否執(zhí)行網(wǎng)格重新劃分。導(dǎo)葉運動過程中,網(wǎng)格質(zhì)量判斷準則為同時滿足網(wǎng)格單元質(zhì)量大于0.3且網(wǎng)格偏斜角小于85°。若不滿足該條件,則首先獲取當(dāng)前時間步的導(dǎo)葉位置,然后經(jīng)過導(dǎo)葉與背景區(qū)域兩者之間的布爾運算獲得新的導(dǎo)葉位置,最后以相同的網(wǎng)格設(shè)置參數(shù)進行網(wǎng)格重新劃分。

    圖4 導(dǎo)葉運動控制流程圖

    圖5 活動導(dǎo)葉運動示意圖

    2.3 求解設(shè)置本文基于ANSYS Fluent進行水輪機負荷變化的數(shù)值計算,湍流模型采用SSTk-ω模型[28-29],該模型改進了近壁區(qū)低雷諾數(shù)流動的處理方式,能較好處理帶有逆壓梯度、流動分離的復(fù)雜流動問題。壁面邊界層則采用全y+壁面函數(shù)處理方式,這種方式的優(yōu)點在于使用混合壁面函數(shù),對精細網(wǎng)格仿效低y+壁面處理,而對粗糙網(wǎng)格仿效高y+壁面處理。通過剛體運動來處理轉(zhuǎn)輪域網(wǎng)格的旋轉(zhuǎn)運動,所有壁面均設(shè)置為無滑移,對流項采用二階迎風(fēng)格式,時間離散則采用高精度的二階格式,收斂標準設(shè)置為殘差小于0.0001。按照M?ssinger等[30]和Sotoudeh等[31]的研究,針對Francis99模型瞬態(tài)數(shù)值研究,瞬態(tài)時間步長設(shè)置為小于轉(zhuǎn)輪旋轉(zhuǎn)4°所用的時間能夠獲得較為可靠的數(shù)值結(jié)果。為了進一步提高數(shù)值精度,捕捉負荷變化過程中的壓力脈動及渦流特性,本文采用對應(yīng)于轉(zhuǎn)輪旋轉(zhuǎn)2°所用的時間,經(jīng)計算瞬態(tài)時間步長為0.001 s,且每個時間步內(nèi)迭代15次。首先進行BEP工況下的穩(wěn)態(tài)計算,并以此為初場進行BEP工況下的非穩(wěn)態(tài)計算,以轉(zhuǎn)輪旋轉(zhuǎn)4圈的非穩(wěn)態(tài)計算結(jié)果作為負荷減小瞬態(tài)計算的初始值。數(shù)值計算進口給定總壓,出口指定靜壓,這種邊界條件的設(shè)置有效保證了數(shù)值計算與模型試驗水輪機進出口壓差的一致性。按照NTNU提供的進口靜壓P1、流量Q和凈水頭H的試驗數(shù)據(jù),可通過式(1)至(3)獲得蝸殼進口總壓Pin和尾水管出口靜壓條件Pout,

    (1)

    (2)

    (3)

    式中:P1為蝸殼進口靜壓,Pa;Pout為尾水管出口靜壓,Pa;Q為流量,m3/s;ρ為密度,取999.8 kg/m3;g為重力加速度,取9.82 m/s2;Ain為蝸殼進口橫截面面積;Aout為尾水管出口橫截面面積;vin為蝸殼進口速度,m/s;vout為尾水管出口速度,m/s;z為水輪機進出口高程差,m。

    此外,本文對NTNU提供的進出口壓力數(shù)據(jù)進行非線性曲線擬合處理,取其擬合結(jié)果作為瞬態(tài)模擬的進出口邊界條件,如圖6中點劃線所示。

    圖6 進出口壓力邊界條件

    2.4 計算域離散及網(wǎng)格無關(guān)性驗證本文采用多面體網(wǎng)格對計算域進行離散,多面體網(wǎng)格的最大優(yōu)點在于其一個控制體周圍有大約10個相鄰單元,因而能以相對較少的控制體獲得更加精確的梯度和局部流動分離解[32-33]。圖7所示為計算域網(wǎng)格示意圖,其中固壁面如導(dǎo)葉、轉(zhuǎn)輪葉片等通過棱柱層網(wǎng)格進行邊界加密處理。

    圖7 計算域多面體網(wǎng)格示意圖

    為消除網(wǎng)格數(shù)量對計算結(jié)果的影響,本文采用美國機械工程師協(xié)會(ASME)推薦的網(wǎng)格收斂指數(shù)(GCI)進行網(wǎng)格離散誤差的估計[34-36]。GCI是一個具有95%置信區(qū)間、表示兩個對比網(wǎng)格中更密網(wǎng)格與漸進值之間距離的指標,并且預(yù)測進一步的網(wǎng)格細化對求解的影響。GCI網(wǎng)格無關(guān)性驗證需要3套不同數(shù)目的網(wǎng)格,分別為細密網(wǎng)格(Fine)、中等網(wǎng)格(Medium)和粗糙網(wǎng)格(Coarse),計算的近似相對誤差ea,外推相對誤差eext以及網(wǎng)格收斂指數(shù)GCI:

    (4)

    (5)

    (6)

    式中:φ為計算所選擇的關(guān)鍵變量;r為網(wǎng)格加密因子;p為采用定點迭代法計算的表觀級數(shù);下標1、2、3分別對應(yīng)于網(wǎng)格方案Fine、Medium、Coarse;上標21和32分別表示網(wǎng)格Fine相對于Medium以及網(wǎng)格Medium相對于Coarse的計算值。

    表1列出了BEP工況下數(shù)值計算離散誤差的計算過程及結(jié)果統(tǒng)計。其中,N1、N2、N3為三種不同密度下的網(wǎng)格數(shù),分別對應(yīng)Fine、Medium和Coarse網(wǎng)格方案。本文選擇水輪機水頭H(單位為m)和轉(zhuǎn)輪扭矩T(單位為N·m)作為網(wǎng)格無關(guān)性測試的關(guān)鍵變量。由表1統(tǒng)計結(jié)果知,三種密度網(wǎng)格的數(shù)值解以漸進形式收斂,表明網(wǎng)格加密有利于平均流場的求解。對Fine和Medium網(wǎng)格而言,計算的水頭不確定度分別為0.072%和0.18%,轉(zhuǎn)輪扭矩不確定度為0.066%和0.25%,表明Fine和Medium網(wǎng)格方案獲得的數(shù)值解誤差較小。為了平衡計算精度與計算資源之間的關(guān)系,本文最終選擇了Medium網(wǎng)格進行數(shù)值計算,網(wǎng)格總數(shù)目為1451萬,其中蝸殼和固定導(dǎo)葉網(wǎng)格數(shù)為292萬,活動導(dǎo)葉為284萬,轉(zhuǎn)輪為514萬,尾水管為361萬。圖8顯示了Medium網(wǎng)格方案BEP工況轉(zhuǎn)輪、尾水管的y+值分布云圖,由圖中可知,轉(zhuǎn)輪壁面最大y+值為11.893,尾水管壁面最大y+值為7.061,且轉(zhuǎn)輪壁面y+值主要范圍在2.5以下,而尾水管壁面y+值主要范圍在1.1左右,可以認為基本適用于SSTk-ω湍流模型。

    表1 數(shù)值計算離散誤差及不確定性統(tǒng)計表

    圖8 轉(zhuǎn)輪及尾水管壁面y+分布

    3 結(jié)果與討論

    3.1 外特性分析為驗證本文所采用數(shù)值策略的可靠性,圖9顯示了數(shù)值模擬與試驗測試獲得的水輪機外特性水頭H和流量Q的對比結(jié)果。

    圖9 負荷減小過程水輪機水頭與流量變化

    圖9結(jié)果顯示,數(shù)值預(yù)測的水輪機水頭與試驗測試趨勢平均值一致性較高。負荷減小過程中,水頭在BEP工況下保持恒定,導(dǎo)葉閉合后,水頭快速上升。t=3.5 s導(dǎo)葉停止閉合,水輪機開始運行在PL工況,此時水頭繼續(xù)上升,在t=4.2 s達到最高值,隨后緩慢下降。水頭變化結(jié)果表明,導(dǎo)葉的閉合對水輪機水頭有一定的提升,且導(dǎo)葉停止運動后這種影響仍會在一定時間內(nèi)持續(xù)。數(shù)值預(yù)測的流量與試驗測試結(jié)果之間存在一定差異,數(shù)值結(jié)果顯示為與導(dǎo)葉開度變化規(guī)律一致的線性變化,而試驗值則為非線性變化規(guī)律,且導(dǎo)葉閉合過程t=1.0 s至t=3.5 s時段內(nèi),試驗記錄的流量值變化率小于數(shù)值結(jié)果。造成試驗與數(shù)值之間差異的主要原因在于試驗測試的流量計存在一定的時間延遲[27]。盡管數(shù)值解與試驗測試存在一定差異,但可通過對比負荷減小過程中定導(dǎo)葉開度BEP階段、PL階段流量值與對應(yīng)穩(wěn)態(tài)試驗測試結(jié)果之間的差異來量化本文數(shù)值流量結(jié)果的可信度。負荷減小過程中,BEP與PL工況數(shù)值計算獲得的流量分別為0.193 24 m3/s和0.136 27 m3/s,對應(yīng)穩(wěn)態(tài)試驗條件流量分別為0.199 59 m3/s和0.139 62 m3/s,僅僅相差3.2%和2.4%,表明數(shù)值計算的流量準確性較高。此外,水輪機流量調(diào)節(jié)是通過改變活動導(dǎo)葉開度實現(xiàn)的,流量變化趨勢與導(dǎo)葉開度密切相關(guān),本文中數(shù)值模擬與試驗測試均顯示了隨導(dǎo)葉閉合流量逐漸減小的變化趨勢,因此可認為本文數(shù)值計算的流量變化規(guī)律具有較高的可信度。綜上所述,采用本文數(shù)值計算的模擬策略,能較好地反映負荷減小過程中水輪機外特性參數(shù)的動態(tài)響應(yīng)。

    3.2 壓力脈動分析為進一步分析負荷減小過程中水輪機內(nèi)部壓力脈動特性,圖10顯示了無葉區(qū)和尾水管內(nèi)VL02、DT5與DT6的壓力脈動系數(shù)Cp,其中Cp為按照式(7)計算的無量綱數(shù),

    (7)

    分析圖10結(jié)果可知,數(shù)值模擬獲得的無葉區(qū)壓力脈動變化規(guī)律與試驗測試結(jié)果吻合地較好,壓力系數(shù)在定導(dǎo)葉開度BEP及PL工況下顯示為均勻波動,而隨著導(dǎo)葉的閉合壓力系數(shù)線性減小。幅值方面,壓力系數(shù)數(shù)值結(jié)果相對于試驗值偏低,這與試驗測試過程中的機械振動、電磁振蕩以及數(shù)值模擬采用的網(wǎng)格、湍流模型和時間步長等有關(guān)。此外,由初始至終止時刻,數(shù)值模擬預(yù)測的壓力系數(shù)幅值保持一致,而試驗結(jié)果在負荷減小過程及PL工況下的壓力脈動幅值要更高。

    圖10 數(shù)值與試驗壓力系數(shù)對比

    尾水管內(nèi)的壓力系數(shù)結(jié)果顯示,數(shù)值模擬低估了尾水管錐管段的壓力脈動幅值,然而,其壓力系數(shù)變化趨勢與試驗值一致性較高,因而在定性分析尾水管內(nèi)的流動特性上具有較高的可靠性。BEP工況下,轉(zhuǎn)輪出口接近法向出流,數(shù)值壓力系數(shù)與試驗壓力系數(shù)平均值均保持為一個固定值,導(dǎo)葉閉合后,壓力系數(shù)迅速下降,隨后隨導(dǎo)葉開度的減小逐漸增加。與BEP工況相似,導(dǎo)葉停止運動進入PL工況,壓力系數(shù)首先迅速上升,隨后緩慢下降,數(shù)值結(jié)果在時刻t=7.2 s以后出現(xiàn)一個顯著的低頻、高幅值周期性脈動,而相同的周期性脈動在試驗過程中約在t=6.5 s后出現(xiàn)。低頻、高幅值周期性壓力系數(shù)的出現(xiàn),與部分負荷工況下尾水管內(nèi)螺旋形尾水管渦帶的出現(xiàn)有關(guān)。

    3.3 轉(zhuǎn)輪軸向水推力分析軸向水推力不僅是水輪機組推力軸承載荷的重要組成部分,而且其脈動特性對會直接引起廠房結(jié)構(gòu)的振動,因而是水輪機設(shè)計和運行過程中的一個重要參數(shù)[37]。為進一步分析軸向水推力相對水輪機負荷變化的響應(yīng)關(guān)系,圖11顯示了負荷減小過程中水輪機軸向水推力隨時間變化規(guī)律,其中正值表示數(shù)值模擬獲得的軸向水推力與參考坐標系z軸方向一致。圖中結(jié)果顯示,BEP工況下軸向水推力絕對值最高且保持在一個常數(shù)附近。軸向水推力變化與導(dǎo)葉運動規(guī)律比較一致,即隨著導(dǎo)葉的閉合而減小,且在導(dǎo)葉開始運動后,軸向水推力快速降低,而在導(dǎo)葉停止閉合后進入PL工況,軸向水推力首先快速上升,然后在一定時間范圍內(nèi)緩慢下降,最后做周期性脈動。

    圖11 負荷減小過程中軸向水推力隨時間變化曲線

    為確定影響軸向水推力的重要因素,本文通過相干分析方法建立軸向水推力與不穩(wěn)定壓力脈動之間的聯(lián)系[38],如圖12所示,圖中分別顯示了軸向水推力與無葉區(qū)測點VL02、轉(zhuǎn)輪域測點P42以及尾水管測點DT6的相干系數(shù)。圖12中橫坐標f/fn定義為一個無量綱頻率比,為計算頻率f與水輪機轉(zhuǎn)頻fn之間的比值。縱坐標為相干系數(shù),用來衡量兩個變量之間的相關(guān)程度,相干系數(shù)越大,表明兩個變量之間關(guān)聯(lián)程度越高。分析圖12可知,除較低頻率范圍及導(dǎo)葉、葉片通過頻率28fn和30fn附近,軸向力相對于DT6的相干性系數(shù)均高于VL02和P42。相干系數(shù)大意味著兩個信號之間的關(guān)聯(lián)程度高,因此認為軸向力Fz與尾水管中DT6的壓力脈動的相關(guān)度要顯著高于無葉區(qū)和轉(zhuǎn)輪內(nèi),表明尾水管內(nèi)部非穩(wěn)態(tài)壓力波動向上游的傳播是造成轉(zhuǎn)輪軸向力波動的重要因素。

    圖12 軸向力與測點VL02、P42、DT6壓力信號相關(guān)性

    3.4 尾水管內(nèi)速度變化情況分析螺旋形尾水管渦帶的形成,與錐管段內(nèi)的軸向速度及切向速度密切相關(guān)[39-40]。圖13和圖14分別為瞬態(tài)過程尾水管L1至L3上的軸向和切向速度分布,圖中軸向速度沿z軸正向為負,切向速度沿轉(zhuǎn)輪轉(zhuǎn)向為正,圖中橫坐標r/r0定義為一個無量綱距離比,為速度測量線L1、L2、L3上坐標位置相對于該線段中間位置距離的比值。由于NTNU未提供切向速度數(shù)據(jù),因此本文僅給出數(shù)值結(jié)果。為便于對不同工況的速度變量進行分析,本文將速度變化分為3個階段,其中BEP工況為階段Ⅰ,導(dǎo)葉開度變化過程BEP至PL為階段Ⅱ,PL工況為階段Ⅲ。

    圖13結(jié)果顯示,水平線L1及L2上數(shù)值預(yù)測的軸向速度與試驗結(jié)果比較吻合,豎直線L3上數(shù)值與試驗結(jié)果存在一定差異,但二者之間的變化趨勢比較一致。第Ⅰ階段,轉(zhuǎn)輪出口水流接近法向出流,軸向速度與主流方向相同,泄水錐下方L1與L2上的軸向速度高于尾水管錐管段壁面附近,而L3沿流向保持一致。導(dǎo)葉閉合的第Ⅱ階段,L1與L2上測得的軸向速度在-0.5

    圖13 軸向速度分布

    分析圖14結(jié)果可知,L1與L2上的切向速度變化比較一致,在第Ⅰ階段均與轉(zhuǎn)輪轉(zhuǎn)向一致,且轉(zhuǎn)軸中心處絕對值較小。第Ⅱ階段,轉(zhuǎn)軸中心附近發(fā)生切向速度的反轉(zhuǎn),且隨時間推移反向切向速度范圍逐漸擴大。第Ⅲ階段,轉(zhuǎn)軸中心附近的正向切向速度區(qū)隨導(dǎo)葉的閉合逐漸增加,而在這個正向速度區(qū)外側(cè),均顯示為反向切向速度。之后,切向速度也表現(xiàn)以轉(zhuǎn)軸為中心線,正向和負向交替出現(xiàn)的切向速度區(qū)。豎直線L3上的切向速度在第I和第Ⅱ階段變化較小,而在第Ⅲ階段后期切向速度出現(xiàn)方向的交替變化,且靠近轉(zhuǎn)輪中心測點經(jīng)歷更多的切向速度方向變化。

    圖14 切向速度分布

    3.5 尾水管渦帶形成分析由BEP工況過渡至PL工況,水輪機尾水管內(nèi)的壓力、速度均顯示了一定的低頻振蕩,這對尾水管渦帶的形成有一定影響。為進一步深入理解水輪機負荷變化過程中尾水管內(nèi)渦結(jié)構(gòu)演化及尾水管渦帶的形成過程,圖15至圖17分別顯示了階段Ⅰ至階段Ⅲ尾水管內(nèi)軸向速度云圖、速度矢量分布及以Q準則顯示的渦結(jié)構(gòu),圖中Q準則閾值q=700 s-2,軸向速度為負時,表示軸向速度方向與參考坐標系z軸相反,與主流方向一致。由于BEP工況不同時刻尾水管內(nèi)流場結(jié)構(gòu)比較接近,因此本文只給出了t=1.0 s時刻的結(jié)果。

    圖15 階段Ⅰ尾水管軸向速度及渦結(jié)構(gòu)

    圖16 階段Ⅱ尾水管軸向速度及渦結(jié)構(gòu)

    圖17 階段Ⅲ尾水管軸向速度及渦結(jié)構(gòu)

    圖15速度分布結(jié)果顯示,BEP工況下尾水管內(nèi)為幾乎與主流方向相同的純軸向流動的穩(wěn)定流場,軸向速度在錐管段中心附近較高,而在錐管段壁面附近較低。由于泄水錐的壁面效應(yīng),在其正下方出現(xiàn)范圍較小的局部回流區(qū)??梢暬瘻u結(jié)構(gòu)等值面結(jié)果未出現(xiàn)自由渦結(jié)構(gòu),僅僅顯示為附著于尾水管進口、泄水錐表面以及彎肘段的片狀渦結(jié)構(gòu),這主要是由于固壁面附近速度梯度較高所致,而且這種片狀渦結(jié)構(gòu)在整個模擬過程的不同時刻均存在。因此,BEP工況下尾水管內(nèi)流動光順,未誘發(fā)不穩(wěn)定渦結(jié)構(gòu)。

    第Ⅱ階段,水輪機活動導(dǎo)葉由9.84°線性關(guān)閉至6.72°,水輪機流量逐漸減小,轉(zhuǎn)輪進口沖角增大。分析圖16結(jié)果可知,t=1.8 s,尾水管內(nèi)的軸向速度由于流量的降低而減小,特別是轉(zhuǎn)軸附近區(qū)域,此時水流方向仍然為軸向流動,未出現(xiàn)回流及分離區(qū)。然而,尾水管錐管段內(nèi)出現(xiàn)了軸對稱的自由渦結(jié)構(gòu),表明水輪機負荷的減小,激發(fā)了尾水管內(nèi)部的不穩(wěn)定流動。t=2.7 s,隨著流量的減小,尾水管內(nèi)的軸向速度進一步降低,錐管段出口與肘管段進口中心區(qū)域出現(xiàn)由肘管段外側(cè)指向內(nèi)側(cè)的流動偏移區(qū)。此時,尾水管內(nèi)的自由渦結(jié)構(gòu)強度顯著增加,由泄水錐延伸至肘管段進口處。導(dǎo)葉繼續(xù)閉合至t=3.0 s,流動偏移強度迅速增強,其中軸向速度流線方向與尾水管軸線呈較大角度,其范圍擴大至彎肘段中部。

    此外,泄水錐下游出現(xiàn)反向軸向速度而形成回流區(qū),此時的渦強度相對t=2.7 s有所增強,但仍然保持為對稱結(jié)構(gòu)。導(dǎo)葉終止閉合時刻t=3.5 s,尾水管內(nèi)的軸向速度持續(xù)降低,流動偏移區(qū)域與回流區(qū)域范圍均顯著擴大,表明水輪機流量的減小使尾水管內(nèi)的水力不穩(wěn)定性進一步增強。然而,圖16中顯示的渦結(jié)構(gòu)沿軸向迅速收縮,而在徑向其渦核半徑增加。綜上所述,水輪機負荷變化過程中,尾水管內(nèi)的軸向速度逐漸降低,誘發(fā)不穩(wěn)定的流動偏移及形成回流區(qū)。可視化渦結(jié)構(gòu)首先顯示為強度較低的軸對稱結(jié)構(gòu),隨后渦強度沿軸向及徑向均顯著增強,最后沿軸向迅速收縮而沿徑向強度增加,而在負荷變化的整個過程中渦形態(tài)均顯示為對稱結(jié)構(gòu)。

    第Ⅲ階段,水輪機運行在定導(dǎo)葉開度PL工況下。t=3.7 s,尾水管內(nèi)的軸向速度分布與t=3.5 s時刻比較接近,但轉(zhuǎn)軸中心區(qū)域的回流強度提高且軸對稱渦結(jié)構(gòu)沿軸向收縮。t=5.5 s時,尾水管錐管段及彎肘段的低速區(qū)范圍顯著增加,對應(yīng)地,尾水管內(nèi)的流動偏移區(qū)域與回流區(qū)域范圍被擴大,且在轉(zhuǎn)軸中心左側(cè)形成兩個明顯的漩渦。然而,相對于t=3.7 s,該時刻的自由渦強度被顯著抑制,僅僅殘留附著在泄水錐表面的部分渦結(jié)構(gòu)。t=7.0 s,尾水管錐管段低速區(qū)內(nèi)出現(xiàn)了多個漩渦區(qū),渦結(jié)構(gòu)強度較上一時刻顯示了較大的差別,該時刻的渦強度有所增強,但渦形態(tài)演化為雙螺旋狀偏心結(jié)構(gòu)。t=8.5 s,錐管段流態(tài)紊亂程度加劇,轉(zhuǎn)軸中心的低速區(qū)不再關(guān)于尾水管中心線對稱,其在尾水管內(nèi)顯示為“S”狀,且靠近壁面的低速區(qū)處出現(xiàn)明顯的漩渦結(jié)構(gòu),對水流形成一定的擠壓作用而使壁面附近的速度提高。Q準則顯示的渦形態(tài)為強度較高、單一的螺旋狀偏心渦結(jié)構(gòu),為典型的部分負荷尾水管渦帶,t=8.5 s之后,該螺旋狀偏心渦帶做與轉(zhuǎn)輪轉(zhuǎn)向相同的旋轉(zhuǎn)運動且其形態(tài)保持固定。綜合以上分析可知,PL工況下,尾水管內(nèi)的回流范圍顯著增加,渦形態(tài)首先為軸對稱形,然后沿軸向劇烈收縮,隨后單一的柱狀渦結(jié)構(gòu)演化為強度升高的雙螺旋狀渦結(jié)構(gòu),最終該雙螺旋狀渦結(jié)構(gòu)合并為單一的螺旋狀偏心渦帶并保持其形態(tài)。

    4 結(jié)論

    本文采用基于多面體網(wǎng)格的動網(wǎng)格技術(shù)對Francis99高水頭混流式模型水輪機負荷減小過渡過程進行了數(shù)值計算,分析了負荷變化對水輪機外特性參數(shù)、不同位置壓力脈動及尾水管內(nèi)渦流結(jié)構(gòu)的影響,主要結(jié)論如下:

    (1)基于多面體網(wǎng)格技術(shù)和進出口壓力邊界條件的水輪機瞬態(tài)數(shù)值研究,能夠較可靠地模擬負荷變化過程中水輪機水頭、流量的變化規(guī)律。負荷變化過程中,無葉區(qū)壓力信號變化趨勢與導(dǎo)葉閉合規(guī)律一致。導(dǎo)葉閉合及停止運動的瞬間對尾水管壓力信號有較大的擾動作用。導(dǎo)葉停止運動進入部分負荷工況,尾水管內(nèi)逐漸產(chǎn)生低頻周期性壓力信號,且轉(zhuǎn)輪軸向力與尾水管內(nèi)壓力脈動關(guān)聯(lián)度最高,表明尾水管內(nèi)的壓力脈動向上游的傳播是影響軸向水推力形成的重要因素。

    (2)導(dǎo)葉閉合導(dǎo)致水輪機流量的降低,對尾水管內(nèi)的橫向流動及回流有一定的促進作用,加劇了尾水管內(nèi)的不穩(wěn)定流動。最優(yōu)工況下,水輪機內(nèi)部為光順流動,尾水管內(nèi)無自由渦結(jié)構(gòu)。負荷減小過程中,尾水管內(nèi)形成的強度較小的軸對稱渦帶結(jié)構(gòu),首先沿軸向拉伸,隨后沿軸向劇烈收縮但徑向尺寸增加。部分負荷工況下,軸對稱渦帶首先潰滅收縮,隨后演變?yōu)楦街谛顾F上的雙螺旋渦帶結(jié)構(gòu),最后雙螺旋渦帶合并為單一的強度較高的螺旋形尾水管渦帶。

    猜你喜歡
    導(dǎo)葉轉(zhuǎn)輪水輪機
    水輪機過流面非金屬材料的修復(fù)及防護
    大中型水斗式水輪機的關(guān)鍵技術(shù)
    詞語大轉(zhuǎn)輪
    ——“AABC”和“無X無X”式詞語
    寺廟里有座大書架——神奇的轉(zhuǎn)輪藏
    水輪機虛擬仿真動畫制作的研究
    壓氣機導(dǎo)葉用GTD222合金鑄造技術(shù)研究
    模具制造(2019年4期)2019-06-24 03:36:50
    我國第一臺分半鑄造的30萬千瓦水輪發(fā)電機轉(zhuǎn)輪實制成功
    蘭臺世界(2017年12期)2017-06-22 14:17:25
    水輪機過流部件改造與節(jié)能增效
    混流式模型水輪機導(dǎo)葉水力矩試驗
    一種新型的電子式導(dǎo)葉主令控制器的設(shè)計及應(yīng)用
    91国产中文字幕| 亚洲av片天天在线观看| 中文字幕制服av| 亚洲成人国产一区在线观看| 精品福利永久在线观看| 成人三级做爰电影| 欧美国产精品va在线观看不卡| ponron亚洲| 亚洲成av片中文字幕在线观看| 中国美女看黄片| 色综合欧美亚洲国产小说| 美女福利国产在线| 每晚都被弄得嗷嗷叫到高潮| 欧美激情极品国产一区二区三区| 桃红色精品国产亚洲av| 法律面前人人平等表现在哪些方面| 香蕉久久夜色| 久久久久久人人人人人| 国产精品一区二区在线观看99| 激情在线观看视频在线高清 | 天天添夜夜摸| 亚洲av第一区精品v没综合| 老司机午夜十八禁免费视频| 欧美乱码精品一区二区三区| 伊人久久大香线蕉亚洲五| 午夜福利乱码中文字幕| 国产精品久久久人人做人人爽| 国产在线一区二区三区精| 精品国产亚洲在线| 亚洲欧美精品综合一区二区三区| 亚洲 国产 在线| 成年人免费黄色播放视频| 十八禁人妻一区二区| 欧美另类亚洲清纯唯美| 色老头精品视频在线观看| 男女午夜视频在线观看| 久久99一区二区三区| 少妇的丰满在线观看| 国产亚洲精品一区二区www | 欧美国产精品一级二级三级| 欧美日韩精品网址| 久久精品国产亚洲av香蕉五月 | 身体一侧抽搐| 岛国毛片在线播放| 日本一区二区免费在线视频| 天天躁日日躁夜夜躁夜夜| 丝瓜视频免费看黄片| 美女高潮到喷水免费观看| 久久久国产成人免费| 国产在线精品亚洲第一网站| 少妇粗大呻吟视频| 91精品三级在线观看| 91国产中文字幕| 精品一区二区三区av网在线观看| 国产一区二区三区视频了| 美女高潮喷水抽搐中文字幕| 久99久视频精品免费| 国产激情欧美一区二区| 日韩免费av在线播放| 在线天堂中文资源库| 久久久国产成人免费| 亚洲 欧美一区二区三区| 少妇粗大呻吟视频| 丝袜美腿诱惑在线| 最新在线观看一区二区三区| 久久精品国产综合久久久| 国产成人精品久久二区二区免费| 女性生殖器流出的白浆| 免费观看人在逋| 久久人妻福利社区极品人妻图片| 成在线人永久免费视频| 亚洲一卡2卡3卡4卡5卡精品中文| 国产成人av激情在线播放| 麻豆成人av在线观看| 国产高清激情床上av| 久久精品国产亚洲av高清一级| 国产国语露脸激情在线看| 美女扒开内裤让男人捅视频| 搡老熟女国产l中国老女人| 少妇 在线观看| 日本五十路高清| 美女高潮喷水抽搐中文字幕| 国产精品久久视频播放| 首页视频小说图片口味搜索| 国产日韩一区二区三区精品不卡| 国产亚洲欧美精品永久| 人人澡人人妻人| 欧美乱码精品一区二区三区| 老熟妇乱子伦视频在线观看| 丰满饥渴人妻一区二区三| 国产无遮挡羞羞视频在线观看| 免费黄频网站在线观看国产| 人妻一区二区av| 狠狠狠狠99中文字幕| 国产精品成人在线| 18禁裸乳无遮挡动漫免费视频| 91精品国产国语对白视频| 亚洲欧美精品综合一区二区三区| 国产男靠女视频免费网站| www.999成人在线观看| 亚洲一码二码三码区别大吗| 很黄的视频免费| 怎么达到女性高潮| 亚洲色图av天堂| 伊人久久大香线蕉亚洲五| 亚洲av美国av| 在线看a的网站| 99香蕉大伊视频| 国产精品久久电影中文字幕 | 国产精品久久久久久精品古装| 精品人妻在线不人妻| 在线看a的网站| 久久久久精品国产欧美久久久| 亚洲成人免费电影在线观看| 亚洲av电影在线进入| 在线av久久热| 成年人午夜在线观看视频| 另类亚洲欧美激情| 在线观看免费午夜福利视频| 免费日韩欧美在线观看| 在线永久观看黄色视频| 欧美日韩视频精品一区| 日韩免费高清中文字幕av| 欧美在线黄色| 青草久久国产| 美国免费a级毛片| 精品国产一区二区三区四区第35| 波多野结衣一区麻豆| 亚洲中文日韩欧美视频| 国产高清激情床上av| 在线观看免费日韩欧美大片| 亚洲一区二区三区不卡视频| 国产精品免费一区二区三区在线 | 一级黄色大片毛片| 国产一区二区三区在线臀色熟女 | 久久精品亚洲熟妇少妇任你| 国产精品二区激情视频| 精品久久蜜臀av无| 身体一侧抽搐| 亚洲一码二码三码区别大吗| 国产激情欧美一区二区| 久久精品aⅴ一区二区三区四区| 18禁裸乳无遮挡动漫免费视频| 美国免费a级毛片| 女人被躁到高潮嗷嗷叫费观| 一区福利在线观看| 亚洲视频免费观看视频| 亚洲九九香蕉| 桃红色精品国产亚洲av| 女人被狂操c到高潮| 亚洲午夜理论影院| 黄色毛片三级朝国网站| 脱女人内裤的视频| 精品国产乱子伦一区二区三区| 嫁个100分男人电影在线观看| 国产亚洲欧美精品永久| 两个人看的免费小视频| 成年人午夜在线观看视频| 老熟女久久久| 岛国毛片在线播放| 国产色视频综合| 欧美丝袜亚洲另类 | 久久久久国产一级毛片高清牌| 黑人猛操日本美女一级片| 在线观看免费高清a一片| 久久香蕉精品热| 精品国产一区二区三区四区第35| 精品久久久久久久毛片微露脸| 成人手机av| 国产精品 欧美亚洲| 久久精品亚洲精品国产色婷小说| 99久久人妻综合| 国产三级黄色录像| 黄片大片在线免费观看| 免费久久久久久久精品成人欧美视频| 咕卡用的链子| 日韩欧美国产一区二区入口| 大陆偷拍与自拍| 大片电影免费在线观看免费| 男人舔女人的私密视频| 中国美女看黄片| 成人精品一区二区免费| 在线永久观看黄色视频| 首页视频小说图片口味搜索| 美女视频免费永久观看网站| 老鸭窝网址在线观看| 精品少妇久久久久久888优播| 亚洲成人免费av在线播放| 最近最新免费中文字幕在线| 在线观看舔阴道视频| 精品国产美女av久久久久小说| 久99久视频精品免费| 免费在线观看视频国产中文字幕亚洲| 亚洲美女黄片视频| 在线看a的网站| 大陆偷拍与自拍| netflix在线观看网站| 久久婷婷成人综合色麻豆| 国产午夜精品久久久久久| 精品无人区乱码1区二区| 男男h啪啪无遮挡| 十八禁高潮呻吟视频| 一个人免费在线观看的高清视频| 最新在线观看一区二区三区| 超碰成人久久| 国产免费av片在线观看野外av| 欧美精品人与动牲交sv欧美| 欧美大码av| 国产精品久久久久久精品古装| 老司机影院毛片| 美女 人体艺术 gogo| 欧美国产精品va在线观看不卡| 久久亚洲精品不卡| 免费看十八禁软件| 人人澡人人妻人| 日韩免费av在线播放| 欧美黑人欧美精品刺激| 国产欧美日韩一区二区精品| 国产精品电影一区二区三区 | 日韩欧美免费精品| 国产成人一区二区三区免费视频网站| 国产单亲对白刺激| 色综合婷婷激情| 亚洲欧美激情在线| 精品一品国产午夜福利视频| 一级,二级,三级黄色视频| 大片电影免费在线观看免费| 在线十欧美十亚洲十日本专区| 国产精品影院久久| 日韩视频一区二区在线观看| 久久这里只有精品19| 亚洲av成人一区二区三| 三上悠亚av全集在线观看| 久久久久久久午夜电影 | e午夜精品久久久久久久| 美国免费a级毛片| 美女午夜性视频免费| 免费人成视频x8x8入口观看| 久久久国产一区二区| 亚洲国产精品合色在线| 久热这里只有精品99| 久久久久久久久久久久大奶| 国产aⅴ精品一区二区三区波| 国产无遮挡羞羞视频在线观看| 亚洲成av片中文字幕在线观看| 亚洲人成电影免费在线| 成人影院久久| 一边摸一边抽搐一进一小说 | 又紧又爽又黄一区二区| 亚洲片人在线观看| 久久亚洲真实| 国产aⅴ精品一区二区三区波| 婷婷精品国产亚洲av在线 | 国产精品免费视频内射| 在线永久观看黄色视频| 精品免费久久久久久久清纯 | 国产精品久久久av美女十八| 男人的好看免费观看在线视频 | 一夜夜www| 久久久久国产一级毛片高清牌| 亚洲一码二码三码区别大吗| 夜夜躁狠狠躁天天躁| 午夜福利视频在线观看免费| 大香蕉久久成人网| netflix在线观看网站| 精品一品国产午夜福利视频| 99热国产这里只有精品6| 欧美激情久久久久久爽电影 | 中出人妻视频一区二区| 欧美人与性动交α欧美软件| 国产又色又爽无遮挡免费看| 欧美日韩视频精品一区| 欧美 亚洲 国产 日韩一| 成年动漫av网址| 后天国语完整版免费观看| 精品卡一卡二卡四卡免费| 日本欧美视频一区| tocl精华| www日本在线高清视频| 两个人免费观看高清视频| 一级作爱视频免费观看| 女性被躁到高潮视频| 十八禁网站免费在线| 成熟少妇高潮喷水视频| ponron亚洲| 亚洲性夜色夜夜综合| 欧美成人免费av一区二区三区 | 日韩熟女老妇一区二区性免费视频| 热99re8久久精品国产| 欧美国产精品va在线观看不卡| 亚洲久久久国产精品| 精品欧美一区二区三区在线| 丰满饥渴人妻一区二区三| 男女床上黄色一级片免费看| 久久久国产成人精品二区 | 国产在线一区二区三区精| 亚洲国产精品一区二区三区在线| aaaaa片日本免费| 一区二区三区国产精品乱码| 国产不卡av网站在线观看| 中文字幕人妻熟女乱码| 色老头精品视频在线观看| 成人亚洲精品一区在线观看| 国产精品永久免费网站| 精品亚洲成a人片在线观看| 视频在线观看一区二区三区| 免费观看精品视频网站| videos熟女内射| 欧美日韩视频精品一区| 亚洲专区字幕在线| videos熟女内射| 久久亚洲真实| 看片在线看免费视频| 国产不卡av网站在线观看| 欧美丝袜亚洲另类 | 午夜福利欧美成人| 国产成人影院久久av| 欧美最黄视频在线播放免费 | 日韩免费高清中文字幕av| 国产精品香港三级国产av潘金莲| 中亚洲国语对白在线视频| 日韩欧美在线二视频 | 制服诱惑二区| 亚洲国产欧美网| 亚洲精品久久成人aⅴ小说| 久9热在线精品视频| 亚洲专区字幕在线| 亚洲视频免费观看视频| 人人妻人人澡人人看| 免费在线观看视频国产中文字幕亚洲| 一进一出好大好爽视频| 日本一区二区免费在线视频| 又黄又爽又免费观看的视频| 亚洲男人天堂网一区| aaaaa片日本免费| 九色亚洲精品在线播放| 高清视频免费观看一区二区| 女人精品久久久久毛片| 国产亚洲精品久久久久久毛片 | 国产日韩一区二区三区精品不卡| 热99久久久久精品小说推荐| 高清在线国产一区| 精品国产一区二区三区四区第35| 熟女少妇亚洲综合色aaa.| 免费黄频网站在线观看国产| 大码成人一级视频| 亚洲国产毛片av蜜桃av| 国产亚洲欧美98| 韩国精品一区二区三区| 99久久精品国产亚洲精品| av线在线观看网站| 国产97色在线日韩免费| 日本欧美视频一区| 可以免费在线观看a视频的电影网站| 午夜免费成人在线视频| 亚洲免费av在线视频| 国产精品影院久久| 在线av久久热| 99re在线观看精品视频| 亚洲综合色网址| 国产真人三级小视频在线观看| 中出人妻视频一区二区| 亚洲片人在线观看| 久久ye,这里只有精品| 精品人妻在线不人妻| 法律面前人人平等表现在哪些方面| 精品国产一区二区三区四区第35| 色综合欧美亚洲国产小说| 在线视频色国产色| 大码成人一级视频| 午夜福利欧美成人| av一本久久久久| 国产精品国产av在线观看| 欧美日韩亚洲国产一区二区在线观看 | 亚洲精品国产区一区二| 久久午夜亚洲精品久久| 婷婷成人精品国产| 国产极品粉嫩免费观看在线| 亚洲精品在线观看二区| 国产精品亚洲av一区麻豆| 午夜福利视频在线观看免费| 亚洲第一欧美日韩一区二区三区| 欧美 日韩 精品 国产| 99国产精品一区二区蜜桃av | 在线观看免费午夜福利视频| netflix在线观看网站| 亚洲美女黄片视频| 欧美日韩视频精品一区| 欧美乱色亚洲激情| 国产91精品成人一区二区三区| 免费不卡黄色视频| 黄网站色视频无遮挡免费观看| 丝袜美腿诱惑在线| 国产日韩一区二区三区精品不卡| 男女午夜视频在线观看| 波多野结衣一区麻豆| 欧美不卡视频在线免费观看 | 午夜免费成人在线视频| 嫁个100分男人电影在线观看| 精品卡一卡二卡四卡免费| 99热只有精品国产| 中文字幕制服av| 久99久视频精品免费| 丝瓜视频免费看黄片| aaaaa片日本免费| 免费黄频网站在线观看国产| 欧美中文综合在线视频| 少妇被粗大的猛进出69影院| 黄色怎么调成土黄色| 美女高潮到喷水免费观看| 亚洲国产欧美网| 免费观看人在逋| 亚洲午夜精品一区,二区,三区| 亚洲黑人精品在线| 在线观看免费视频日本深夜| 欧美日韩视频精品一区| 99国产极品粉嫩在线观看| 首页视频小说图片口味搜索| 国产高清激情床上av| 制服人妻中文乱码| 久久香蕉国产精品| 咕卡用的链子| 欧美精品人与动牲交sv欧美| 建设人人有责人人尽责人人享有的| 成人黄色视频免费在线看| 国产欧美日韩一区二区三区在线| 欧美丝袜亚洲另类 | 国产亚洲精品第一综合不卡| 国产亚洲欧美精品永久| 欧美黄色片欧美黄色片| 热re99久久国产66热| 国产亚洲精品久久久久5区| 极品少妇高潮喷水抽搐| 国产麻豆69| 久久久国产欧美日韩av| 国产精品免费大片| 亚洲av美国av| 精品卡一卡二卡四卡免费| 天堂中文最新版在线下载| 国产精品一区二区精品视频观看| 亚洲色图 男人天堂 中文字幕| 一进一出抽搐gif免费好疼 | 超碰成人久久| 丝袜人妻中文字幕| 91字幕亚洲| 亚洲一区高清亚洲精品| 丰满的人妻完整版| 日日摸夜夜添夜夜添小说| 精品国产一区二区三区久久久樱花| 亚洲va日本ⅴa欧美va伊人久久| 久久久久久久精品吃奶| 妹子高潮喷水视频| 日本wwww免费看| 欧美在线黄色| 久久热在线av| 国产成人精品在线电影| 欧美激情久久久久久爽电影 | 母亲3免费完整高清在线观看| 国产亚洲av高清不卡| 精品熟女少妇八av免费久了| 国产在线观看jvid| 国产高清国产精品国产三级| 国产一区二区激情短视频| 成人永久免费在线观看视频| 在线观看午夜福利视频| 伊人久久大香线蕉亚洲五| 亚洲精品乱久久久久久| 亚洲第一青青草原| 午夜日韩欧美国产| 老鸭窝网址在线观看| 欧美亚洲日本最大视频资源| 国产区一区二久久| 久久精品国产清高在天天线| 久久 成人 亚洲| 51午夜福利影视在线观看| 一边摸一边做爽爽视频免费| 日本wwww免费看| 国产不卡一卡二| 精品久久久久久久久久免费视频 | 色精品久久人妻99蜜桃| 纯流量卡能插随身wifi吗| 首页视频小说图片口味搜索| 18禁国产床啪视频网站| 99久久99久久久精品蜜桃| 满18在线观看网站| 国产99白浆流出| 亚洲成a人片在线一区二区| 一区二区日韩欧美中文字幕| 丁香六月欧美| 亚洲一码二码三码区别大吗| 久久婷婷成人综合色麻豆| 亚洲免费av在线视频| 他把我摸到了高潮在线观看| 成人精品一区二区免费| 99国产精品一区二区蜜桃av | 亚洲精品乱久久久久久| 精品乱码久久久久久99久播| 咕卡用的链子| 日韩大码丰满熟妇| 老熟女久久久| 高清在线国产一区| 人人妻人人澡人人爽人人夜夜| 色老头精品视频在线观看| 后天国语完整版免费观看| 久久国产乱子伦精品免费另类| 国产片内射在线| 国产免费现黄频在线看| 国产精品1区2区在线观看. | 身体一侧抽搐| 久久草成人影院| 一本一本久久a久久精品综合妖精| 亚洲欧美一区二区三区久久| 一级,二级,三级黄色视频| 男人舔女人的私密视频| 国产成人啪精品午夜网站| 一二三四社区在线视频社区8| 高清欧美精品videossex| 在线观看免费日韩欧美大片| 无遮挡黄片免费观看| 制服人妻中文乱码| 亚洲九九香蕉| 热re99久久精品国产66热6| 亚洲一区高清亚洲精品| 日韩欧美三级三区| 午夜福利影视在线免费观看| 国产精品成人在线| 亚洲国产精品sss在线观看 | 两性午夜刺激爽爽歪歪视频在线观看 | 国产成人av教育| 国产高清激情床上av| 99国产精品99久久久久| 电影成人av| 自线自在国产av| 欧美国产精品一级二级三级| 国产精品一区二区在线观看99| 亚洲专区国产一区二区| 亚洲欧美激情综合另类| 久久中文看片网| 国产成+人综合+亚洲专区| 青草久久国产| 国产又爽黄色视频| 窝窝影院91人妻| 啪啪无遮挡十八禁网站| 丝瓜视频免费看黄片| 高清在线国产一区| 国产av精品麻豆| 国产亚洲一区二区精品| 成年女人毛片免费观看观看9 | 久久久久精品国产欧美久久久| 午夜福利免费观看在线| 侵犯人妻中文字幕一二三四区| 亚洲熟女精品中文字幕| 精品一区二区三卡| 91老司机精品| 男女免费视频国产| 国产精品综合久久久久久久免费 | 又大又爽又粗| 精品亚洲成国产av| 黄色怎么调成土黄色| 99国产精品一区二区蜜桃av | 免费久久久久久久精品成人欧美视频| 日韩熟女老妇一区二区性免费视频| 亚洲欧美精品综合一区二区三区| aaaaa片日本免费| 免费在线观看完整版高清| 无限看片的www在线观看| 亚洲熟女精品中文字幕| 在线观看免费视频网站a站| 看黄色毛片网站| 1024视频免费在线观看| 久久久国产欧美日韩av| 脱女人内裤的视频| 国产精品 国内视频| 十八禁高潮呻吟视频| 丰满迷人的少妇在线观看| 日本a在线网址| 精品电影一区二区在线| 欧美精品av麻豆av| 免费女性裸体啪啪无遮挡网站| 不卡一级毛片| 校园春色视频在线观看| 国产精品.久久久| 国产视频一区二区在线看| 久久精品国产99精品国产亚洲性色 | 王馨瑶露胸无遮挡在线观看| 在线观看免费日韩欧美大片| 亚洲国产看品久久| 亚洲av成人一区二区三| 久久午夜亚洲精品久久| 欧美黑人精品巨大| 女人爽到高潮嗷嗷叫在线视频| 男人操女人黄网站| 成年动漫av网址| 咕卡用的链子| 69精品国产乱码久久久| av视频免费观看在线观看| 成人永久免费在线观看视频| 99香蕉大伊视频| 国产不卡av网站在线观看| 久久精品亚洲熟妇少妇任你| 国产欧美日韩综合在线一区二区| 精品国产美女av久久久久小说| 国产精品久久久久成人av| 国产精品美女特级片免费视频播放器 | 中文字幕最新亚洲高清| 一夜夜www| 欧美+亚洲+日韩+国产| 制服人妻中文乱码| 最近最新中文字幕大全免费视频| 日韩欧美国产一区二区入口| 丰满人妻熟妇乱又伦精品不卡| 久久人妻福利社区极品人妻图片| 亚洲国产欧美日韩在线播放| 丰满人妻熟妇乱又伦精品不卡| 免费女性裸体啪啪无遮挡网站|