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

    10 MW固態(tài)燃料釷基熔鹽堆穩(wěn)態(tài)物理-熱工耦合

    2017-08-25 08:58:08陳家豪1張海青1朱智勇1
    核技術(shù) 2017年8期
    關(guān)鍵詞:功率密度熔鹽熱工

    陳家豪1,2 張海青1 朱智勇1

    ?

    10 MW固態(tài)燃料釷基熔鹽堆穩(wěn)態(tài)物理-熱工耦合

    陳家豪張海青朱智勇

    1(中國科學(xué)院上海應(yīng)用物理研究所嘉定園區(qū) 上海 201800) 2(中國科學(xué)院大學(xué) 北京 100049)

    固態(tài)燃料釷基熔鹽堆(Thorium Molten Salt Reactor-Solid Fuel, TMSR-SF1)作為第四代先進(jìn)核反應(yīng)堆堆型之一,繼承了熔鹽冷卻劑和球形燃料元件的許多優(yōu)點(diǎn)和技術(shù)基礎(chǔ),具有良好的經(jīng)濟(jì)性、設(shè)計(jì)上的固有安全性、釷鈾燃料的可持續(xù)性和防核擴(kuò)散性。本文以10 MW固態(tài)燃料釷基熔鹽堆為模型,利用MCNP (Monte Carlo N Particle Transport Code)和ANSYS Fluent等模擬程序?qū)ζ溥M(jìn)行多物理耦合分析,同時(shí)利用C++語言編寫了堆芯活性區(qū)的物理-熱工耦合計(jì)算程序,實(shí)現(xiàn)了MCNP計(jì)算結(jié)果與Fluent程序的對(duì)接,并且通過對(duì)比耦合前后結(jié)果,分析了堆芯功率密度分布、有效增殖因子、溫度分布等主要參數(shù),為熔鹽堆的設(shè)計(jì)、安全性評(píng)估和操作運(yùn)行提供了參考依據(jù)。

    固態(tài)燃料球,釷基熔鹽堆,物理熱工耦合,穩(wěn)態(tài)分析

    固態(tài)熔鹽堆是由美國科學(xué)家于2003年提出的第四代反應(yīng)堆堆型之一,其主要特點(diǎn)包括:1) 使用氟化物熔鹽冷卻(熔鹽堆);2) 使用三向同性(Tri-structural Isotropic, TRISO)包覆顆粒燃料(高溫氣冷堆)。包覆顆粒燃料元件有球形和柱形兩種基本類型,其特點(diǎn)是能夠在高溫下運(yùn)行,具有本征的安全性。固態(tài)燃料釷基熔鹽堆(Thorium Molten Salt Reactor-Solid Fuel, TMSR-SF1)是中國科學(xué)院“未來先進(jìn)核裂變能”戰(zhàn)略性先導(dǎo)科技專項(xiàng)優(yōu)先選擇并設(shè)計(jì)的兩種反應(yīng)堆堆型之一。2012年4月,中國科學(xué)院TMSR專項(xiàng)啟動(dòng)第一個(gè)固態(tài)燃料釷基熔鹽堆的設(shè)計(jì)和建造。該堆型采用球形燃料元件,利用熔鹽冷卻劑把固態(tài)燃料球產(chǎn)生的裂變能有效傳遞到熱力系統(tǒng)進(jìn)行能量轉(zhuǎn)換。

    核反應(yīng)堆的運(yùn)行是多物理場(chǎng)相互作用的過程,對(duì)核反應(yīng)堆的分析研究已從單純的中子物理學(xué)分析、熱工水力分析、燃料性能分析、結(jié)構(gòu)材料應(yīng)力分析和系統(tǒng)響應(yīng)分析,發(fā)展到如今的多物理場(chǎng)耦合分析。中子物理與熱工耦合計(jì)算在反應(yīng)堆分析中具有重要意義,其結(jié)果可用于更好地理解反應(yīng)堆中的復(fù)雜現(xiàn)象。耦合計(jì)算時(shí),分別使用中子物理計(jì)算程序和熱工水力計(jì)算程序計(jì)算相同的模型問題并在迭代過程中交換數(shù)據(jù),直到結(jié)果收斂達(dá)到穩(wěn)定狀態(tài)。

    針對(duì)多種熔鹽堆堆型的物理分析已有較多研究。Cheng等針對(duì)熔鹽快堆開發(fā)了三維穩(wěn)態(tài)模擬程序,并用于模擬熔鹽快堆的穩(wěn)態(tài)運(yùn)行狀況。薛春等則設(shè)計(jì)了一種將燃料球裝入燃料組件的熔鹽堆堆芯并對(duì)其進(jìn)行了設(shè)計(jì)計(jì)算。Kophazi等通過耦合 DALTON 和 THERM 程序,構(gòu)建了一個(gè)三維耦合程序,并應(yīng)用于 MSRE 穩(wěn)態(tài)和瞬態(tài)多物理耦合分析。Zhang等則針對(duì)2 MW熔鹽實(shí)驗(yàn)堆的屏蔽層計(jì)算了中子注量分布和溫度分布。Guo等采用 MCNP (Monte Carlo N Particle Transport Code)與多通道熱工水力程序耦合的方式分析了簡化的MSRE堆芯穩(wěn)態(tài)功率分布、溫度分布、壓降和流量分配,并對(duì)堆芯設(shè)計(jì)進(jìn)行了優(yōu)化。Zhou等采用MCNP與CFD (Computational Fluid Dynamics)熱工水力程序耦合的方法,模擬了幾種堆芯條件驟變時(shí)液態(tài)燃料熔鹽堆的瞬態(tài)情況。而Li等通過將PB-AHTR (Pebble Bed Advanced High Temperature Reactor)堆芯模型作為基礎(chǔ),利用清華大學(xué)自主開發(fā)的RMC (Reactor Monte Carlo codes)中子輸運(yùn)程序和熱工水力程序CFX耦合,分析了固態(tài)燃料釷基熔鹽堆穩(wěn)態(tài)條件下該方法的準(zhǔn)確性及耦合效率。何杰等采用MCNP和ANSYS Fluent兩種程序,利用Python語言編寫了后處理程序MTF (MCNP to Fluent)并對(duì)液態(tài)熔鹽堆TMSR-LF1的主屏蔽體開展了耦合計(jì)算分析。本文以TMSR-SF1為主要研究對(duì)象,以中子輸運(yùn)程序MCNP和熱工水力程序ANSYS Fluent分別開展了TMSR-SF1的三維中子物理和三維熱工水力模型計(jì)算,通過C++編程語言進(jìn)行上述兩個(gè)模擬程序的數(shù)據(jù)交換以達(dá)到多物理耦合的目的,并對(duì)模擬結(jié)果中耦合前后的部分中子物理和熱工水力指標(biāo)進(jìn)行了分析比較。本文對(duì)球床堆堆芯活性區(qū)的TRISO顆粒模型進(jìn)行了簡化以提高計(jì)算效率,并利用燃料球的堆積模塊對(duì)球床進(jìn)行了合理的劃分,使球床堆的耦合數(shù)據(jù)交換得以進(jìn)行。

    1 計(jì)算方法及模型

    1.1 堆芯結(jié)構(gòu)

    TMSR-SF1堆芯包括堆芯活性區(qū)和石墨反射層,如圖1所示,其中活性區(qū)由燃料球隨機(jī)堆積的球床和流經(jīng)球床的熔鹽組成,熔鹽冷卻劑在燃料球間的間隙中流動(dòng)。上下圓臺(tái)和側(cè)面的環(huán)繞石墨構(gòu)建了反射層,中間圓柱區(qū)為活性區(qū),在反射層設(shè)有中子源、停堆系統(tǒng)、實(shí)驗(yàn)測(cè)量等專用通道。燃料元件為直徑6 cm的固態(tài)燃料球,其中燃料UO中U富集度為17%,每個(gè)燃料球裝有7.0 g鈾。一回路冷卻劑采用含高富集度Li (99.995%)的FLiBe (2LiF-BeF)熔鹽。

    以圖1堆芯模型為基礎(chǔ),建立中子物理模型用于蒙特卡羅方法計(jì)算,其中包括堆芯活性區(qū)、上下及側(cè)面反射層、冷卻劑腔室、控制棒通道等。

    圖1 TMSR-SF1堆芯示意圖 (a) 橫截面,(b) 豎截面

    1.2 堆芯物理計(jì)算

    由于TMSR-SF1所裝載的一萬多個(gè)燃料球由基體石墨和包含有數(shù)千個(gè)TRISO燃料顆粒組成,而TRISO顆粒由4層包覆層和燃料核心構(gòu)成,對(duì)如此眾多的燃料球建模并得到精確計(jì)算結(jié)果是非常耗時(shí)的。因此,F(xiàn)ratoni等提出通過將TRISO顆粒中的包覆層與基體石墨按照原子數(shù)之比合算成一種材料并保留燃料核心,從而達(dá)到簡化模型結(jié)構(gòu)、縮短計(jì)算時(shí)間的目的。表1是對(duì)TMSR-SF1所使用燃料球的計(jì)算結(jié)果對(duì)比,從表1中發(fā)現(xiàn),該方法可用于TMSR-SF1并有效地縮短計(jì)算時(shí)間。

    表1 同質(zhì)化TRISO顆粒效果

    由于蒙特卡羅方法在計(jì)算有效增殖因子和反應(yīng)率計(jì)數(shù)時(shí)存在概率性誤差,因此在進(jìn)行耦合模擬前,對(duì)該模型進(jìn)行了所需中子數(shù)及計(jì)算誤差分析。本文計(jì)算中每次迭代中子數(shù)為50 000,迭代次數(shù)為550次,跳過前50次迭代,計(jì)算的標(biāo)準(zhǔn)誤差為0.000 16,統(tǒng)計(jì)誤差小于0.1,計(jì)算時(shí)間約4 h。

    進(jìn)行中子物理計(jì)算所用截面數(shù)據(jù)庫是在ENDF/B-VII數(shù)據(jù)庫基礎(chǔ)上,采用了針對(duì)熔鹽高溫特點(diǎn)加工的數(shù)據(jù)庫ACEDATA,此截面數(shù)據(jù)庫選取了300?2200K之間共30個(gè)溫度點(diǎn),通過NJOY加工而成。堆芯活性區(qū)軸向向上劃分為11層,徑向按照燃料球體心立方堆積模塊,每層劃分為88份,總共836個(gè)柵元,并對(duì)其使用F7計(jì)數(shù)卡計(jì)數(shù)以供后續(xù)功率密度計(jì)算。

    1.3 熱工水力計(jì)算

    使用CFD方法進(jìn)行針對(duì)球床堆堆芯的熱工水力模擬研究中,球床的復(fù)雜性經(jīng)常造成所需網(wǎng)格數(shù)龐大且難以得到精準(zhǔn)的計(jì)算結(jié)果。Wu等針對(duì)高溫氣冷堆HTGR (High Temperature Gas-cooled Reactor)的球床模型,運(yùn)用CFD軟件對(duì)比了真實(shí)燃料球間隙通道和多孔介質(zhì)模型通道的計(jì)算結(jié)果,得出了多孔介質(zhì)模型可用于計(jì)算堆芯活性區(qū)中冷卻劑整體溫度、流速、沿程壓降、密度等主要參數(shù)的結(jié)論。因此,本文采用多孔介質(zhì)模型描述燃料球球床結(jié)構(gòu),以簡化網(wǎng)格、節(jié)省計(jì)算資源。

    1.3.1 模型介紹

    熱工水力模型根據(jù)中子物理模型的堆芯活性區(qū)冷卻劑流動(dòng)區(qū)域建立,包括上下石墨反射層的冷卻劑通道及堆芯活性區(qū),具體結(jié)構(gòu)如圖2所示。該模型使用ANSYS ICEM劃分網(wǎng)格,并與物理模型的相關(guān)區(qū)域在空間上完全重合。

    圖2 熱工水力模型示意圖

    1.3.2 邊界條件

    1) 下反射層冷卻劑入口通道邊界條件。冷卻劑進(jìn)入下反射層入口通道,由下至上經(jīng)過堆芯活性區(qū)并通過上反射層冷卻劑通道離開堆芯。本文將各冷卻劑入口通道流速設(shè)置為TMSR-SF1概念設(shè)計(jì)報(bào)告中參數(shù),平均流速0.2235 m?s,最大流速0.744m?s,最小流速0.033 m?s,冷卻劑入口溫度設(shè)置為873.15 K。

    2) 內(nèi)熱源加載。堆芯活性區(qū)內(nèi)熱源通過將MCNP計(jì)算得到的沉積能量轉(zhuǎn)換為對(duì)應(yīng)空間的功率密度,再通過耦合程序?qū)氲皆O(shè)定熱工水力模型功率密度、空間坐標(biāo)信息的UDF (User-defined Function)中,內(nèi)熱源設(shè)置根據(jù)不同坐標(biāo)位置加載對(duì)應(yīng)的功率密度。

    3) 其他邊界條件。出口邊界條件設(shè)置為壓力出口??紤]重力影響,設(shè)置重力加速度為9.8 m?s。采用SIMPLE算法二階迎風(fēng)格式。

    1.4 耦合過程

    冷卻劑密度隨溫度的變化、核反應(yīng)截面隨溫度的變化以及堆芯結(jié)構(gòu)材料的熱膨脹是中子物理模擬中最為重要的溫度反饋信息。本文耦合過程中熱工水力模擬僅考慮溫度及冷卻劑密度變化的影響。

    1.4.1 數(shù)據(jù)交換

    通過外耦合的方式耦合MCNP和ANSYS Fluent。在每一耦合循環(huán)中,兩個(gè)模型分別進(jìn)行計(jì)算,并通過C++語言開發(fā)程序完成數(shù)據(jù)交換過程。堆芯活性區(qū)在兩個(gè)模型中皆在軸向分隔為11層;在每一次耦合循環(huán)中,耦合程序?qū)⒅凶游锢砟P椭?1層的沉積能量數(shù)據(jù)轉(zhuǎn)換為功率密度并用于更新熱工水力模型中內(nèi)熱源設(shè)置的更新。耦合步驟概括為如下步驟,圖3顯示了模擬數(shù)據(jù)在兩程序之間的傳遞過程。

    圖3 MCNP與ANSYS Fluent耦合流程圖

    1) MCNP以冷態(tài)堆芯為初始條件進(jìn)行中子物理計(jì)算,從輸出文件提取F7卡對(duì)各柵元的平均裂變沉積能計(jì)數(shù),通過式(1)轉(zhuǎn)為功率密度:

    式中:P是第個(gè)柵元的功率密度,W·m;(F)是第個(gè)柵元裂變沉積能量,MeV;V是第個(gè)柵元的體積,m;是反應(yīng)堆堆芯熱功率,W。

    2) 以此功率密度分布更新ANSYS Fluent內(nèi)熱源設(shè)置,并進(jìn)行熱工水力計(jì)算,得到堆芯活性區(qū)溫度分布和密度分布。

    3) 以溫度和密度分布數(shù)據(jù)更新MCNP輸入文件中的柵元密度信息和調(diào)取的截面庫,再次計(jì)算。

    4) 比較前后兩次中子物理計(jì)算結(jié)果中的有效增殖因子大小,判斷耦合過程是否收斂。當(dāng)偏差較大時(shí),重復(fù)2)?4)步驟的操作。

    1.4.2 核反應(yīng)截面的多普勒展寬

    多物理耦合的重要環(huán)節(jié)之一即更新中子物理計(jì)算的核反應(yīng)截面。本文所采用的截面庫更新方法為偽材料法。通過偽材料法更新截面庫,可大幅減少計(jì)算所需截面庫的存儲(chǔ)體積,且無需預(yù)處理截面庫,節(jié)省計(jì)算資源;此方法的缺點(diǎn)是在更新MCNP計(jì)算的截面庫時(shí),并不是直接產(chǎn)生目標(biāo)溫度點(diǎn)所對(duì)應(yīng)的截面庫,而是通過計(jì)算目標(biāo)溫度點(diǎn)所在溫度區(qū)間的上下兩組截面庫的比例分?jǐn)?shù),近似擬合截面庫數(shù)據(jù),且無法擬合熱中子散射截面庫(,)數(shù)據(jù),需要選擇與目標(biāo)溫度最接近的溫度點(diǎn)的熱中子散射截面。

    各核素的核反應(yīng)截面自ACEDATA提取,在MCNP計(jì)算溫度信息更新后,由式(2)計(jì)算新的溫度點(diǎn)所需提取截面庫的分?jǐn)?shù):

    式中:為新的柵元溫度;和分別為截面庫中所在區(qū)間的最高溫度和最低溫度;為最低溫度對(duì)應(yīng)的核截面所占份額。所求柵元溫度對(duì)應(yīng)的截 面為:

    (3)

    2 計(jì)算結(jié)果及討論

    2.1 計(jì)算驗(yàn)證

    本文所使用的中子物理模型計(jì)算每次迭代使用50000個(gè)粒子,計(jì)算結(jié)果中有效增殖因子為1.04031,概念設(shè)計(jì)報(bào)告計(jì)算結(jié)果為1.03990,計(jì)算誤差為0.0394%。

    熱工計(jì)算結(jié)果精度雖足夠高,但在實(shí)際工程應(yīng)用中達(dá)到毫米量級(jí)的流速變化對(duì)主要參數(shù)無明顯影響,因此保留三位有效數(shù)字。熱工水力模型耦合前后入口截面平均速度均為0.233 m?s,耦合情況下出口截面平均流速為0.215 m?s。未耦合情況下出口截面平均溫度為903.935 K,耦合情況下出口截面平均溫度為904.035 K,計(jì)算誤差為0.01%。

    耦合程序的驗(yàn)證采用Hu等針對(duì)使用MCNP和Fluent進(jìn)行中子物理與熱工水力耦合可行性研究所采用的模型,對(duì)比手動(dòng)輸入功率密度與耦合程序轉(zhuǎn)換的計(jì)算結(jié)果來驗(yàn)證計(jì)算方法是否正確。

    驗(yàn)證模型由10 cm×10 cm×10 cm的正方體按4×4×4的方式堆積而成,共64塊。左側(cè)32個(gè)方塊為U富集度3%的UO,右側(cè)為輕水,MCNP模型與此對(duì)應(yīng)。冷卻劑流動(dòng)方向?yàn)?軸方向,入口流速為5 m?s。手動(dòng)輸入與程序輸入的計(jì)算結(jié)果符合較好,驗(yàn)證了耦合程序的準(zhǔn)確性和可靠性。

    2.2 物理耦合計(jì)算結(jié)果分析

    TMSR-SF1堆芯經(jīng)過耦合計(jì)算的有效增值因子如圖4所示。從圖4中可以看出,未耦合情況下為1.04031,經(jīng)過一次耦合循環(huán)后,有效增殖因子相對(duì)于未耦合計(jì)算結(jié)果有明顯變化,增大約1.08%,之后該因子不再隨著耦合循環(huán)次數(shù)增加而變化,在1.051附近波動(dòng);耦合后,有效增殖因子明顯增大的原因是MCNP輸入文件中的材料截面庫根據(jù)溫度反饋更新后,部分位置的實(shí)際溫度與初始冷態(tài)溫度相差較大,材料與粒子相互作用的反應(yīng)截面變化劇烈,最終影響有效增殖因子。而該因子在一定范圍內(nèi)波動(dòng)的情況,則是由蒙特卡羅方法固有的統(tǒng)計(jì)誤差造成的。

    圖4 TMSR-SF1堆芯keff隨耦合次數(shù)的變化

    耦合計(jì)算后TMSR-SF1堆芯活性區(qū)中子能譜如圖5所示,不同軸向高度中子通量密度分布如圖6所示,活性區(qū)功率密度分布如圖7所示。堆芯活性區(qū)平均功率密度5.19 MW?m,未耦合狀態(tài)下最大功率密度9.17 MW?m,耦合狀態(tài)下最大功率密度9.37 MW?m,整個(gè)活性區(qū)耦合前后功率密度改變量為0.24%?3.31%。從圖6、7中可知,耦合前后不同軸向高度的功率密度與中子通量密度變化趨勢(shì)相同,堆芯軸向中心部分的中子通量密度最高,在堆芯活性區(qū)中部最低,靠近反射層處,因部分中子自反射層返回堆芯,使堆芯邊緣功率密度增大。對(duì)比相同軸高耦合前后的徑向功率密度分布發(fā)現(xiàn),功率密度分布一致,但部分位置功率密度值相差較大的原因是未耦合時(shí)MCNP計(jì)算使用平均溫度,耦合后當(dāng)部分位置實(shí)際溫度與平均溫度偏差較大時(shí),功率密度值便產(chǎn)生更大變化。

    圖5 堆芯活性區(qū)中子能譜

    圖6 不同軸向高度中子通量密度分布

    圖7 不同軸向高度功率密度分布

    2.3 熱工水力耦合計(jì)算結(jié)果與分析

    耦合前后堆芯活性區(qū)溫度分布如圖8所示。由圖8可以看出,熔鹽在堆芯活性區(qū)中向上流動(dòng),在球床下方空腔中溫度變化較小,進(jìn)入球床后被燃料球加熱,溫度逐漸升高,在頂部區(qū)域達(dá)到最高,與未耦合情況下堆芯最高溫度912.098 K相比,耦合后變?yōu)?13.793 K,升高1.695K。圖9為出口孔道溫度分布,中心出口孔道溫度較高,外側(cè)孔道溫度較低。

    圖8 活性區(qū)XY截面溫度分布 (a) 未耦合,(b) 耦合

    圖9 堆芯出口截面溫度分布

    采用耦合計(jì)算前后的堆芯徑向溫度場(chǎng)如圖10所示。從圖10中可以看出,未耦合和耦合計(jì)算的徑向溫度在剛進(jìn)入球床(0.55 m)處差別較小,但隨著冷卻劑在堆芯活性區(qū)中向上流動(dòng),溫差逐漸增大,最大溫差出現(xiàn)在堆芯中心軸向高度1.5 m處,為7.584K。兩者的溫度分布變化趨勢(shì)和溫度值有明顯差異,耦合計(jì)算的溫度變化趨勢(shì)與功率密度分布 一致。

    圖10 不同軸向高度的溫度徑向分布

    3 結(jié)語

    本文根據(jù)現(xiàn)階段10 MW TMSR-SF1的主要設(shè)計(jì)參數(shù),建立了該堆型的中子物理及熱工水力模型,并通過C++編程語言實(shí)現(xiàn)了兩個(gè)模型間的耦合模擬,編制了用于在MCNP和ANSYS Fluent之間交換堆芯活性區(qū)功率密度分布、溫度場(chǎng)等數(shù)據(jù)的耦合程序,并分析了穩(wěn)態(tài)情況下的有效增殖因子、堆芯功率密度分布、溫度分布等數(shù)據(jù),主要結(jié)論如下:

    1) TMSR-SF1堆芯在裝料高度為130 cm、控制棒完全抽出的情況下,與未耦合時(shí)單用MCNP的計(jì)算數(shù)值相比,耦合后有效增殖因子從1.0403升至1.0515,增大約1.08%,堆芯功率密度改變范圍為0.24%?3.31%。

    2) 在入口溫度為873.15 K、入口平均流速0.2335m?s的情況下,與未耦合時(shí)單用Fluent的計(jì)算數(shù)值相比,耦合后堆芯最高溫度由912.098 K升至913.793 K,升高1.695K,而耦合前后冷卻劑最大溫差為7.584K,出現(xiàn)在球床中心軸向高度1.5 m處。

    中子物理-熱工耦合研究可提供核反應(yīng)堆運(yùn)行中更為準(zhǔn)確的信息,為研究TMSR-SF1中冷卻劑流動(dòng)狀態(tài)、溫度場(chǎng)、材料性質(zhì)對(duì)反應(yīng)堆運(yùn)行的影響、功率密度分布及核燃料燃耗研究提供了重要參考。

    1 Ingersoll D T, Forsberg C W, Ott L J,. Status of preconceptual design of the advanced high-temperature reactor (AHTR)[M]. United States: Department of Energy, 2004. DOI:10.2172/861752.

    2 江綿恒, 徐洪杰, 戴志敏. 未來先進(jìn)核裂變能——TMSR核能系統(tǒng)[J]. 中國科學(xué)院院刊, 2012, 27(3): 366?374. DOI: 10.3969/j.issn.1000-3045.2012.03.016. JIANG Mianheng, XU Hongjie, DAI Zhimin. Advanced fission energy program - TMSR nuclear energy system[J]. Bulletin of Chinese Academy of Sciences, 2012, 27(3): 366?374. DOI: 10.3969/j.issn.1000-3045.2012.03.016.

    3 Cheng MS,Dai ZM. Development of a three dimension multi-physics code for molten salt fast reactor[J]. Nuclear Science and Techniques, 2014, 25(1): 010601. DOI: 10.13538/j.1001-8042/nst.25.010601.

    4 薛春, 張海青, 朱智勇, 等. 組件型熔鹽堆燃料組件的設(shè)計(jì)研究[J]. 核技術(shù), 2016, 39(9): 090602. DOI: 10.11889/j.02533219.2016.hjs.39.090602.XUE Chun, ZHANG Haiqing, ZHU Zhiyong,. Design of fuel assembly for molten-salt-cooled reactors[J]. Nuclear Techniques, 2016, 39(9): 090602. DOI: 10.11889/ j.0253-3219.2016.hjs.39.090602.

    5 Kophazi J, Lathouwers D, Kloosterman J L. Development of a three-dimensional time-dependent calculation scheme for molten salt reactors and validation of the measurement data of the molten salt reactor experiment[J]. Nuclear Science and Engineering, 2009, 163(2): 118?131. DOI: 10.13182/NSE163-118.

    6 Zhang Z H, Xia X B, Cai J,. Simulation of radiation dose distribution and thermal analysis for the bulk shielding of an optimized molten salt reactor[J]. Nuclear Science and Techniques, 2015, 26(4): 040603. DOI:10.13538/j.1001-8042/nst.26.040603.

    7 Guo Z P, Wang C L, Zhang D L,. The effects of core zoning on optimization of design analysis of molten salt reactor[J]. Nuclear Engineering and Design, 2013, 265: 967?977. DOI: 10.1016/j.nucengdes.2013.09.036.

    8 Guo Z P, Zhou J J, Zhang D L,. Coupled neutronics/thermal-hydraulics for analysis of molten salt reactor[J]. Nuclear Engineering and Design, 2013, 258(2):144?156. DOI:10.1016/j.nucengdes.2013.01.013.

    9 Zhou J J, Wang C L, An H Z,. Three dimensional neutronic/thermal-hydraulic coupled simulation of MSR in steady state condition[J]. Nuclear Engineering and Design, 2014,267(2):88?99.DOI: 10.1016/j.nucengdes.2013.11. 074.

    10 Li L S, Yuan H M, Wang K. Coupling of RMC and CFX for analysis of Pebble Bed-Advanced High Temperature Reactor core[J]. Nuclear Engineering and Design, 2012, 250:385?391. DOI:10.1016/j.nucengdes.2012.05.036.

    11 何杰, 夏曉彬, 蔡軍, 等. 2 MW 液態(tài)釷基熔鹽實(shí)驗(yàn)堆主屏蔽溫度場(chǎng)分析[J]. 核技術(shù), 2016, 39(4): 040601. DOI: 10.11889/j.0253-3219.2016.hjs.39.040601.HE Jie, XIA Xiaobin, CAI Jun,. Temperature field analysis for the main shielding of the 2-MW thorium-based molten salt experimental reactor[J]. Nuclear Techniques, 2016, 39(4): 040601. DOI: 10.11889/j. 0253-3219.2016.hjs.39.040601.

    12 Fratoni M, Greenspan E. Neutronic feasibility assessment of liquid salt-cooled pebble bed reactors[J]. Nuclear Science and Engineering, 2011, 168(1): 1?22. DOI: 10.13182/NSE10-38.

    13 Wu C Y, Ferng Y M, Chieng C C,. Investigating the advantages and disadvantages of realistic approach and porous approach for closely packed pebbles in CFD simulation[J]. Nuclear Engineering and Design, 2010, 240(5): 1151?1159. DOI:10.1016/j.nucengdes.2010.01. 015.

    14 Vazquez M, Tsige-Tamirat H, Ammirabile L,. Coupled neutronics thermal-hydraulics analysis using Monte Carlo and sub-channel codes[J]. Nuclear Engineering and Design, 2012, 250(3):403?411. DOI: 10.1016/j.nucengdes.2012. 06.007.

    15 Conlin J L, Ji W, Lee J C,. Pseudo material construct for coupled neutronic-thermal-hydraulic analysis of VHTGR[J]. Transactions of American Nuclear Society, 2005, 92: 225?227.

    16 Hu J, Rizwan-uddin. Coupled neutronics and thermal-hydraulics simulations using MCNP and FLUENT[J]. Transactions of the American Nuclear Society, 2008, 98:606?608.

    Coupled neutronic and thermal-hydraulic analysis of TMSR-SF1 at steady state

    CHEN JiahaoZHANG HaiqingZHU Zhiyong

    1(Shanghai Institute of Applied Physics, Chinese Academy of Sciences, Jiading Campus, Shanghai 201800, China) 2(University of Chinese Academy of Sciences, Beijing 100049, China)

    Background: Neutronic and thermal-hydraulic simulations of advanced reactors can affect each other’s results. Purpose: This study focuses on coupling neutronic and thermal-hydraulic simulations to achieve more accurate results for future developments of 10-MW solid-fueled thorium molten salt experimental reactor (TMSR-SF1). Methods: A program converting the MCNP (Monte Carlo N particle transport code) results to the spatial distribution of power density within the active region was created using C++ programming language. The spatial distribution data were loaded into the ANSYS Fluent in the form of user-defined function (UDF) to accomplish the coupling of the two simulation processes. In regards of TMSR-SF’s original design parameters, the physical and thermal-hydraulic models of the whole core were established by using MCNP and ANSYS Fluent respectively. Results: The coupling method is feasible and can be used to obtain reliable results. The changes in coolant’s temperature and velocity in the active region are dependent on the power density distribution. The changes in multiplication factor, power density and maximum of discrepancy in coolant temperature are 1.08%, 3.31% and 7.584 K, respectively. Conclusion: It is necessary to take the coupling effects of the reactor core into consideration in the design of associated reactor systems. In addition, the results confirm that the design parameters of the TMSR-SF1 are reasonable.

    Fuel pebble, Thorium molten salt reactor, Neutronics and thermal-hydraulics coupling, Steady-state analysis

    TL426

    10.11889/j.0253-3219.2017.hjs.40.080603

    中國科學(xué)院戰(zhàn)略先導(dǎo)科技專項(xiàng)(No.XDA02030200)資助

    陳家豪,男,1992年出生,2014年畢業(yè)于華中科技大學(xué),現(xiàn)為碩士研究生,研究領(lǐng)域?yàn)楹思夹g(shù)及應(yīng)用

    朱智勇,E-mail: zhuzhiyong@sinap.ac.cn

    2017-03-10,

    2017-06-06

    Strategic Priority Program of Chinese Academy of Sciences (No.XDA02030200)

    CHEN Jiahao, male, born in 1992, graduated from Huazhong University of Science and Technology in 2014, master student, focusing on nuclear technology and application

    ZHU Zhiyong, E-mail: zhuzhiyong@sinap.ac.cn

    2017-03-10, accepted date: 2017-06-06

    猜你喜歡
    功率密度熔鹽熱工
    熔鹽在片堿生產(chǎn)中的應(yīng)用
    NaF-KF熔鹽體系制備Ti2CTx材料的研究
    熱工儀表自動(dòng)化安裝探討的認(rèn)識(shí)
    智能控制在電廠熱工自動(dòng)化中的應(yīng)用
    純鈦的熔鹽滲硼
    智能控制在電廠熱工自動(dòng)化中的應(yīng)用
    大型燃?xì)馊埯}爐的研發(fā)和工藝控制
    工業(yè)爐(2016年1期)2016-02-27 12:34:11
    高效高功率密度低噪聲電機(jī)研究
    PrimePACKTM結(jié)合最新IGBT5和.XT模塊工藝延長產(chǎn)品壽命,提高功率密度
    提高火電廠熱工管理精細(xì)化水平探索
    十分钟在线观看高清视频www| 少妇猛男粗大的猛烈进出视频| 亚洲av成人一区二区三| 亚洲欧美精品综合一区二区三区| 丝瓜视频免费看黄片| 高清视频免费观看一区二区| 成年人黄色毛片网站| 免费av中文字幕在线| 99久久99久久久精品蜜桃| 丝袜在线中文字幕| 一区福利在线观看| 考比视频在线观看| 男女午夜视频在线观看| 中文字幕制服av| 亚洲精品在线美女| 建设人人有责人人尽责人人享有的| 少妇 在线观看| a级毛片黄视频| 一二三四在线观看免费中文在| 欧美精品亚洲一区二区| 亚洲一区中文字幕在线| 9热在线视频观看99| 又紧又爽又黄一区二区| 麻豆乱淫一区二区| 欧美黑人精品巨大| 亚洲av国产av综合av卡| 久久热在线av| 飞空精品影院首页| av有码第一页| 日韩一卡2卡3卡4卡2021年| 亚洲精品av麻豆狂野| 一边摸一边做爽爽视频免费| 99精品欧美一区二区三区四区| 久久精品91无色码中文字幕| 国产精品一区二区免费欧美| 黄网站色视频无遮挡免费观看| 免费黄频网站在线观看国产| 免费高清在线观看日韩| 最新在线观看一区二区三区| 精品国产国语对白av| 亚洲 国产 在线| 一二三四社区在线视频社区8| 高清毛片免费观看视频网站 | 国产xxxxx性猛交| 黄色视频,在线免费观看| 黄色片一级片一级黄色片| 人人澡人人妻人| 亚洲中文字幕日韩| 久久久久精品人妻al黑| 一进一出抽搐动态| 真人做人爱边吃奶动态| 男女午夜视频在线观看| 99精国产麻豆久久婷婷| 一本—道久久a久久精品蜜桃钙片| 亚洲精品av麻豆狂野| 国产欧美日韩一区二区三区在线| 国产成人精品在线电影| 国产精品亚洲av一区麻豆| 欧美乱码精品一区二区三区| 在线观看人妻少妇| 视频区欧美日本亚洲| 一级片免费观看大全| 无人区码免费观看不卡 | 久久精品国产亚洲av香蕉五月 | 9热在线视频观看99| 中文字幕av电影在线播放| 人人澡人人妻人| 丁香六月天网| 久久天堂一区二区三区四区| 久久狼人影院| 日韩视频在线欧美| tube8黄色片| www.自偷自拍.com| 久久狼人影院| 精品免费久久久久久久清纯 | 视频区欧美日本亚洲| 中国美女看黄片| 黄色片一级片一级黄色片| 久久久精品94久久精品| 午夜日韩欧美国产| 一边摸一边抽搐一进一出视频| 妹子高潮喷水视频| 亚洲第一av免费看| 99在线人妻在线中文字幕 | 视频在线观看一区二区三区| 久久这里只有精品19| 丰满迷人的少妇在线观看| 一本综合久久免费| 欧美日韩黄片免| 亚洲美女黄片视频| 亚洲第一青青草原| 免费观看人在逋| 国产av又大| 2018国产大陆天天弄谢| 成人精品一区二区免费| 麻豆成人av在线观看| videos熟女内射| 考比视频在线观看| 丝袜美腿诱惑在线| 亚洲一码二码三码区别大吗| 大香蕉久久成人网| 欧美日本中文国产一区发布| 亚洲av美国av| 亚洲精品美女久久久久99蜜臀| 国产精品久久久人人做人人爽| 欧美黑人欧美精品刺激| 女性生殖器流出的白浆| 一区二区三区精品91| 日韩有码中文字幕| 午夜激情久久久久久久| 国产高清videossex| 亚洲中文日韩欧美视频| 黄色成人免费大全| 嫩草影视91久久| 国产日韩欧美亚洲二区| 亚洲精品美女久久av网站| 一级片免费观看大全| 多毛熟女@视频| 午夜激情av网站| 亚洲视频免费观看视频| 国产精品秋霞免费鲁丝片| 国产高清国产精品国产三级| 黑人巨大精品欧美一区二区mp4| 天天躁狠狠躁夜夜躁狠狠躁| 欧美日韩中文字幕国产精品一区二区三区 | 嫁个100分男人电影在线观看| 久久影院123| 国产高清视频在线播放一区| 精品免费久久久久久久清纯 | 少妇 在线观看| 99re6热这里在线精品视频| 麻豆国产av国片精品| 真人做人爱边吃奶动态| 欧美变态另类bdsm刘玥| 18禁黄网站禁片午夜丰满| 美女福利国产在线| 欧美精品高潮呻吟av久久| 国产成人精品无人区| 国产国语露脸激情在线看| 视频区图区小说| 99精国产麻豆久久婷婷| 一夜夜www| 日本a在线网址| 97在线人人人人妻| 国产野战对白在线观看| 国产不卡av网站在线观看| 亚洲成av片中文字幕在线观看| 日韩欧美免费精品| 欧美乱妇无乱码| 亚洲国产av新网站| 成人18禁高潮啪啪吃奶动态图| 精品国内亚洲2022精品成人 | 国产精品久久电影中文字幕 | 亚洲欧洲精品一区二区精品久久久| 真人做人爱边吃奶动态| 激情在线观看视频在线高清 | 一区二区三区激情视频| 动漫黄色视频在线观看| 亚洲精品在线观看二区| 日韩中文字幕欧美一区二区| 91成人精品电影| 男女午夜视频在线观看| 宅男免费午夜| 黄色视频,在线免费观看| 在线看a的网站| 亚洲黑人精品在线| 18禁裸乳无遮挡动漫免费视频| 精品第一国产精品| 一本综合久久免费| 12—13女人毛片做爰片一| 日韩欧美一区视频在线观看| 日韩人妻精品一区2区三区| 国产亚洲欧美精品永久| 午夜视频精品福利| 男女免费视频国产| 午夜免费鲁丝| 亚洲国产看品久久| 一进一出抽搐动态| 亚洲专区国产一区二区| 麻豆av在线久日| 久久午夜综合久久蜜桃| 国产野战对白在线观看| 国产熟女午夜一区二区三区| 国产精品免费视频内射| 久久人妻福利社区极品人妻图片| 在线播放国产精品三级| 国产男靠女视频免费网站| 久久免费观看电影| 啦啦啦在线免费观看视频4| 天天影视国产精品| 国产精品偷伦视频观看了| 亚洲国产欧美网| 国产高清激情床上av| 精品国产乱码久久久久久小说| 中文字幕精品免费在线观看视频| 亚洲精品粉嫩美女一区| 国产亚洲午夜精品一区二区久久| 99久久人妻综合| 亚洲av欧美aⅴ国产| 99精品在免费线老司机午夜| 一区二区日韩欧美中文字幕| 国产成人av教育| 两性午夜刺激爽爽歪歪视频在线观看 | 午夜91福利影院| 一区二区日韩欧美中文字幕| 成人国产av品久久久| 高清黄色对白视频在线免费看| 黑人欧美特级aaaaaa片| 日韩精品免费视频一区二区三区| 精品久久久久久电影网| 亚洲精品一卡2卡三卡4卡5卡| 久久久水蜜桃国产精品网| 久久久久久亚洲精品国产蜜桃av| 国产精品久久久久久精品古装| 国产亚洲欧美在线一区二区| 成人影院久久| 精品视频人人做人人爽| 久久久久久亚洲精品国产蜜桃av| 老汉色∧v一级毛片| 亚洲av日韩在线播放| 青青草视频在线视频观看| 亚洲国产毛片av蜜桃av| 国产精品成人在线| 手机成人av网站| 国产午夜精品久久久久久| 亚洲精品国产一区二区精华液| 婷婷成人精品国产| 久久久久久久精品吃奶| 咕卡用的链子| 久久久久精品人妻al黑| 久久午夜综合久久蜜桃| 亚洲少妇的诱惑av| 欧美人与性动交α欧美精品济南到| 女性被躁到高潮视频| 国产欧美日韩精品亚洲av| aaaaa片日本免费| 一二三四在线观看免费中文在| 视频区欧美日本亚洲| www.精华液| 91字幕亚洲| 亚洲七黄色美女视频| 精品第一国产精品| 91成人精品电影| 久久久国产成人免费| 久久久久网色| 中文字幕另类日韩欧美亚洲嫩草| 精品高清国产在线一区| 久久久水蜜桃国产精品网| 夫妻午夜视频| 巨乳人妻的诱惑在线观看| 大陆偷拍与自拍| 美女扒开内裤让男人捅视频| 欧美+亚洲+日韩+国产| 亚洲精品国产精品久久久不卡| 成人国语在线视频| 人人妻,人人澡人人爽秒播| 国产精品av久久久久免费| 天天添夜夜摸| 热99国产精品久久久久久7| 欧美av亚洲av综合av国产av| 国产精品免费一区二区三区在线 | 国产淫语在线视频| 狠狠狠狠99中文字幕| 欧美中文综合在线视频| 九色亚洲精品在线播放| 亚洲 国产 在线| 俄罗斯特黄特色一大片| 视频区欧美日本亚洲| 性少妇av在线| 窝窝影院91人妻| 久久亚洲精品不卡| 成年人黄色毛片网站| 黄色毛片三级朝国网站| 99国产极品粉嫩在线观看| 大香蕉久久网| 久久中文字幕一级| 高清av免费在线| 亚洲精品成人av观看孕妇| 无遮挡黄片免费观看| 99国产综合亚洲精品| 国内毛片毛片毛片毛片毛片| 美女主播在线视频| 我的亚洲天堂| 美国免费a级毛片| 男男h啪啪无遮挡| 天天操日日干夜夜撸| 成年人免费黄色播放视频| 久久影院123| 午夜福利影视在线免费观看| 一区二区日韩欧美中文字幕| 免费高清在线观看日韩| 国产极品粉嫩免费观看在线| 多毛熟女@视频| 丁香六月欧美| 搡老熟女国产l中国老女人| 久久久国产欧美日韩av| 日韩三级视频一区二区三区| 国产精品成人在线| 午夜免费鲁丝| 纯流量卡能插随身wifi吗| 19禁男女啪啪无遮挡网站| bbb黄色大片| 亚洲中文日韩欧美视频| 成人av一区二区三区在线看| 人人妻人人爽人人添夜夜欢视频| 国产黄色免费在线视频| 亚洲情色 制服丝袜| 国产精品影院久久| 亚洲av电影在线进入| 18禁黄网站禁片午夜丰满| 我要看黄色一级片免费的| 一个人免费在线观看的高清视频| 一级毛片精品| 亚洲成av片中文字幕在线观看| 国产人伦9x9x在线观看| 亚洲国产看品久久| 精品午夜福利视频在线观看一区 | 高清视频免费观看一区二区| 欧美国产精品一级二级三级| 国产精品香港三级国产av潘金莲| 大码成人一级视频| 97在线人人人人妻| 桃花免费在线播放| av天堂在线播放| 视频区图区小说| 成人av一区二区三区在线看| 99精品久久久久人妻精品| 在线观看www视频免费| 国产精品 国内视频| 免费在线观看视频国产中文字幕亚洲| 极品教师在线免费播放| 美国免费a级毛片| 老司机影院毛片| 国精品久久久久久国模美| 宅男免费午夜| svipshipincom国产片| 日本欧美视频一区| 一本色道久久久久久精品综合| 成人特级黄色片久久久久久久 | 亚洲av片天天在线观看| 久久精品人人爽人人爽视色| 脱女人内裤的视频| 国产三级黄色录像| 免费一级毛片在线播放高清视频 | 9191精品国产免费久久| 一个人免费在线观看的高清视频| 精品国产国语对白av| 成人免费观看视频高清| 啦啦啦中文免费视频观看日本| 亚洲久久久国产精品| 制服人妻中文乱码| 男女床上黄色一级片免费看| 啦啦啦中文免费视频观看日本| 99国产精品免费福利视频| av天堂久久9| 国产一卡二卡三卡精品| 国产aⅴ精品一区二区三区波| 欧美黑人精品巨大| 亚洲欧洲精品一区二区精品久久久| 久久精品国产亚洲av香蕉五月 | 精品少妇一区二区三区视频日本电影| 高清毛片免费观看视频网站 | 91老司机精品| 高清视频免费观看一区二区| 亚洲 欧美一区二区三区| 99在线人妻在线中文字幕 | 欧美 日韩 精品 国产| 日本av免费视频播放| 蜜桃国产av成人99| bbb黄色大片| 亚洲 国产 在线| 人人妻人人澡人人看| 嫩草影视91久久| 精品一区二区三区四区五区乱码| 男女午夜视频在线观看| 久久久久国产一级毛片高清牌| 男女午夜视频在线观看| 精品午夜福利视频在线观看一区 | 日韩一卡2卡3卡4卡2021年| 亚洲伊人色综图| 久久久精品国产亚洲av高清涩受| 精品亚洲成国产av| 亚洲精品乱久久久久久| 在线观看免费视频网站a站| 夫妻午夜视频| 久久久欧美国产精品| 夫妻午夜视频| 18禁国产床啪视频网站| 免费在线观看日本一区| 人人妻人人澡人人爽人人夜夜| 91精品国产国语对白视频| 90打野战视频偷拍视频| 在线观看舔阴道视频| 国产精品av久久久久免费| 无限看片的www在线观看| 国产免费福利视频在线观看| av免费在线观看网站| 免费看十八禁软件| 久久久精品区二区三区| 国产熟女午夜一区二区三区| 乱人伦中国视频| 极品人妻少妇av视频| 91成人精品电影| 亚洲第一欧美日韩一区二区三区 | 成人国产av品久久久| av片东京热男人的天堂| 日韩视频在线欧美| 国产人伦9x9x在线观看| 老汉色av国产亚洲站长工具| 国产成人免费无遮挡视频| 亚洲中文av在线| 黑人操中国人逼视频| 国产亚洲精品久久久久5区| 欧美激情久久久久久爽电影 | 纯流量卡能插随身wifi吗| 热99re8久久精品国产| 新久久久久国产一级毛片| 可以免费在线观看a视频的电影网站| 午夜福利欧美成人| 精品久久蜜臀av无| 亚洲专区中文字幕在线| 国产精品久久久久久人妻精品电影 | 天天影视国产精品| av网站免费在线观看视频| 亚洲情色 制服丝袜| 久久狼人影院| 99热国产这里只有精品6| 日韩中文字幕欧美一区二区| 亚洲性夜色夜夜综合| 国产片内射在线| 人人妻,人人澡人人爽秒播| 亚洲天堂av无毛| 五月天丁香电影| 国产精品九九99| 国产高清国产精品国产三级| 午夜福利欧美成人| 操美女的视频在线观看| 日韩一区二区三区影片| 欧美精品一区二区大全| 咕卡用的链子| 欧美日韩中文字幕国产精品一区二区三区 | 最新美女视频免费是黄的| 久久婷婷成人综合色麻豆| 亚洲一码二码三码区别大吗| 欧美精品av麻豆av| 啦啦啦免费观看视频1| 菩萨蛮人人尽说江南好唐韦庄| 国产精品久久久久久精品古装| 男女边摸边吃奶| av不卡在线播放| 亚洲九九香蕉| 在线观看免费高清a一片| 精品人妻熟女毛片av久久网站| 成人影院久久| 80岁老熟妇乱子伦牲交| 色综合欧美亚洲国产小说| 亚洲美女黄片视频| 菩萨蛮人人尽说江南好唐韦庄| 伊人久久大香线蕉亚洲五| 欧美成狂野欧美在线观看| 97人妻天天添夜夜摸| 黄色视频,在线免费观看| 中文字幕高清在线视频| 日本精品一区二区三区蜜桃| 咕卡用的链子| 亚洲一区二区三区欧美精品| 婷婷成人精品国产| 一本—道久久a久久精品蜜桃钙片| 18禁国产床啪视频网站| 丰满迷人的少妇在线观看| 日韩熟女老妇一区二区性免费视频| 色综合婷婷激情| 国产成人精品久久二区二区91| 亚洲中文字幕日韩| 黄色片一级片一级黄色片| 免费少妇av软件| 最新美女视频免费是黄的| 女人久久www免费人成看片| 久久精品亚洲熟妇少妇任你| 久久性视频一级片| 嫩草影视91久久| av网站在线播放免费| 欧美黄色片欧美黄色片| 男男h啪啪无遮挡| 国产成人精品无人区| 热re99久久精品国产66热6| 极品人妻少妇av视频| 亚洲精品久久成人aⅴ小说| 午夜免费鲁丝| 男女之事视频高清在线观看| 国产亚洲一区二区精品| e午夜精品久久久久久久| 两个人免费观看高清视频| 日日夜夜操网爽| 中文字幕另类日韩欧美亚洲嫩草| 亚洲欧美激情在线| 一进一出抽搐动态| 午夜福利乱码中文字幕| 看免费av毛片| 在线看a的网站| 999久久久精品免费观看国产| 久久久精品国产亚洲av高清涩受| 国产精品影院久久| 免费在线观看黄色视频的| 国产精品1区2区在线观看. | 一区二区三区激情视频| 国产精品 国内视频| 精品福利观看| 亚洲av电影在线进入| 免费一级毛片在线播放高清视频 | 久久人妻福利社区极品人妻图片| 中文字幕人妻丝袜制服| 高潮久久久久久久久久久不卡| 国产av一区二区精品久久| 国产精品 欧美亚洲| 在线观看一区二区三区激情| 青草久久国产| 我的亚洲天堂| 宅男免费午夜| 国产在线精品亚洲第一网站| 欧美黄色片欧美黄色片| 人人澡人人妻人| 十八禁网站网址无遮挡| 啦啦啦中文免费视频观看日本| 亚洲久久久国产精品| 国产99久久九九免费精品| 十八禁网站网址无遮挡| 国产成人免费无遮挡视频| 午夜福利在线免费观看网站| 亚洲欧美一区二区三区久久| 午夜福利在线免费观看网站| 91九色精品人成在线观看| 亚洲九九香蕉| 99国产精品一区二区三区| 亚洲一码二码三码区别大吗| 黄色视频在线播放观看不卡| 日韩人妻精品一区2区三区| 亚洲欧美日韩另类电影网站| 欧美日韩成人在线一区二区| 国产一区二区三区视频了| 女警被强在线播放| 精品视频人人做人人爽| 真人做人爱边吃奶动态| 久久精品91无色码中文字幕| 国产男靠女视频免费网站| 久久久久精品国产欧美久久久| 成人亚洲精品一区在线观看| 国产精品一区二区精品视频观看| 久久这里只有精品19| 一进一出好大好爽视频| 亚洲午夜理论影院| 又大又爽又粗| 黑人巨大精品欧美一区二区蜜桃| 正在播放国产对白刺激| 精品国产乱码久久久久久男人| av免费在线观看网站| 免费人妻精品一区二区三区视频| 欧美亚洲日本最大视频资源| 如日韩欧美国产精品一区二区三区| 久久人人97超碰香蕉20202| 夜夜夜夜夜久久久久| 新久久久久国产一级毛片| av片东京热男人的天堂| 久久久国产欧美日韩av| 久久久国产一区二区| 精品人妻在线不人妻| 最近最新中文字幕大全电影3 | 高清在线国产一区| 日本精品一区二区三区蜜桃| 免费一级毛片在线播放高清视频 | 国产欧美日韩综合在线一区二区| 亚洲第一av免费看| av国产精品久久久久影院| 久久久久久久久久久久大奶| 精品亚洲乱码少妇综合久久| 亚洲精品国产色婷婷电影| 久久人妻福利社区极品人妻图片| 亚洲精品在线美女| 高清av免费在线| 91成年电影在线观看| 九色亚洲精品在线播放| 两人在一起打扑克的视频| 99国产精品免费福利视频| 五月天丁香电影| 一区二区av电影网| 午夜福利,免费看| 亚洲va日本ⅴa欧美va伊人久久| 亚洲人成电影免费在线| 两个人免费观看高清视频| 一本综合久久免费| 亚洲一卡2卡3卡4卡5卡精品中文| 巨乳人妻的诱惑在线观看| 成人精品一区二区免费| 在线观看免费日韩欧美大片| 无遮挡黄片免费观看| 日日爽夜夜爽网站| 丝袜人妻中文字幕| 精品亚洲成a人片在线观看| 国产精品一区二区免费欧美| 亚洲精品中文字幕一二三四区 | 国产激情久久老熟女| 91成年电影在线观看| 妹子高潮喷水视频| 国产亚洲午夜精品一区二区久久| 久久国产精品大桥未久av| 999精品在线视频| 男女高潮啪啪啪动态图| 欧美日本中文国产一区发布| 丰满少妇做爰视频| 久久亚洲真实| 亚洲国产av影院在线观看| a级毛片黄视频| 久久精品91无色码中文字幕| 精品国产超薄肉色丝袜足j| 99热国产这里只有精品6|