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

    First-principles calculations of high pressure and temperature properties of Fe7C3

    2023-09-05 08:48:44LiLiFan范莉莉XunLiu劉勛ChangGao高暢ZhongLiLiu劉中利YanLiLi李艷麗andHaiJunHuang黃海軍
    Chinese Physics B 2023年7期
    關(guān)鍵詞:莉莉海軍

    Li-Li Fan(范莉莉), Xun Liu(劉勛), Chang Gao(高暢), Zhong-Li Liu(劉中利),Yan-Li Li(李艷麗), and Hai-Jun Huang(黃海軍),?

    1School of Science,Wuhan University of Technology,Wuhan 430070,China

    2School of Materials Science and Engineering,Harbin Institute of Technology,Harbin 150001,China

    Keywords: iron carbide,phase stability,thermoelastic properties,sound velocities,inner core

    1.Introduction

    As is well known, the Earth’s core is composed of Fe–Ni alloy with a small quantity of light elements.Based on geochemical and meteorological studies, light elements most likely originate from carbon (C), hydrogen (H), oxygen (O),sulfur (S), and silicon (Si).[1]Until now, the light elemental composition of the Earth’s core has been poorly understood.Among these light elements,C has been suggested as a major light elemental component of the Earth’s core,as C is enriched in the solar system[2]and has a high affinity for liquid iron under core–mantle differentiation conditions.[3]In addition, the shear wave velocity and Poisson’s ratio of Fe7C3are similar to those of the inner core[4,5]and Fe7C3melting experiments indicate that C might be a potential light element in the inner core.[6]

    However, whether Fe7C3is stable under core condition and the nature of its relatively stable structure remain unclear.Based on previous experimental results, two potential stable structures can be reasonably hypothesized at high temperature and pressure: hexagonal structure (space groupP63mc; h-Fe7C3) and orthorhombic structure (space groupPnma,Pmn21mc, orPbca; o-Fe7C3).[5,7–9]Nakajimaet al.found that h-Fe7C3is stable up to 71.5 GPa and 1973 K by usingin situx-ray diffraction.[10]Prescheret al.suggested that o-Fe7C3is the more stable phase up to 205 GPa and 3700 K by using single-crystal and powder x-ray diffraction.[5]Theoretical calculations have been similarly inconclusive and often inconsistent.Fanget al.found that the o-Fe7C3phase is more stable than h-Fe7C3under ambient conditions.[11]Razaet al.found that o-Fe7C3(space groupPbca)is stable below 100 GPa,while h-Fe7C3becomes more stable above 150 GPa and requires higher temperatures for stability.[12]

    The thermal equation of state(EOS)of Fe7C3was measured experimentally only at low pressures and temperatures.Litasovet al.[13]and Nakajimaet al.[10]measured the EOS below 30 GPa; however, their results deviate widely when extrapolated to the high temperature and pressure conditions of the inner core.Chen and Zhang acquired compressional(VP) and shear wave velocity (VS) data for Fe7C3up to 154 GPa at 300 K by using x-ray diffraction.[4]Subsequently,Prescheret al.investigated theVPandVSfor Fe7C3up to 158 GPa and 300 K by using the same method.[5]However,these sound velocity measurements of Fe7C3exhibited large discrepancies.[4,5]Therefore,it is necessary to investigate the crystal structure and thermoelastic properties of Fe7C3under core condition to ascertain the amount of carbon present in the inner core.

    Here, using first-principles calculations, we estimate the free energy,EOS,and sound velocity of Fe7C3.The EOS and sound velocity data for Fe7C3are calculated at 0 K–5000 K,which is close to the inner core boundary temperature obtained by phonon spectrum calculations.We also use our findings to discuss whether carbon is the main light element present in the inner core.

    2.Methods

    All calculations were performed by using the Viennaab initiosimulation package (VASP) based on density functional theory (DFT).[14]The conjugate gradient (CG) algorithm was adopted for the geometric optimization of the crystal structures.Exchange–correlation interactions were considered as the Perdew–Burke–Ernzerhof (PBE) parametrization of the gradient approximation (GGA).[15]The electron–ion interactions were represented by projected augmented-wave(PAW) pseudopotentials.[16,17]The simulation cell contained 40 atoms for o-Fe7C3.Prior to formal calculation, we calculated the free energy values of 40 and 20 atoms for h-Fe7C3at some pressure points, obtaining highly similar results.To balance the computational efficiency and precision, we selected 20 atoms for h-Fe7C3calculation at all pressure points.The cutoff energy was 500 eV.The Monkhorst–Pack mesh was 6×4×7 for o-Fe7C3and 6×6×9 for h-Fe7C3.[18]The phonon spectra and free energy of o-Fe7C3and h-Fe7C3were calculated and analyzed by using the Phonopy[19]and Phasego programs.[20]The phonon spectra and Helmholtz free energy were calculated up to 360 GPa in the quasi-harmonic approximation (QHA).During the thermodynamic simulations, we constructed 1×1×2 and 1×2×2 supercells for o- and h-Fe7C3,and set denserk-point meshes to 7×5×3 and 6×6×9 grids,respectively.[21–23]

    Previous studies have reported that Fe7C3undergoes a magnetic transition from a ferromagnetic (FM) to nonmagnetic (NM) phase at high pressure.[5]Therefore, we performed spin-limited and unlimited calculations to investigate the effect of spin on the EOS.

    The elastic constant of Fe7C3was estimated by using the stress-strain method.[24]Force constant was calculated in real space by using density functional perturbation theory.The bulk modulus and shear modulus were calculated from the elastic tensor by using the Voight–Hill–Reuss method.[24,25]Sound velocities were obtained from the following formulas:

    whereVPis the compressional wave velocity,VSis the shear wave velocity,Bis the bulk modulus,Gis the shear modulus,andρis the density.Considering the symmetry of Fe7C3,BandGcan be expressed as

    whereCx(x=11,22,33,44,55,66,12,13,or 23)values are all elastic constants.

    3.Results

    3.1.Relative stabilities of o-and h-phase Fe7C3

    To determine which structure is stable, we calculate the Gibbs free energy values of o-Fe7C3and h-Fe7C3at high pressure and temperature (Fig.1).The difference in Gibbs free energy in o-Fe7C3and h-Fe7C3is a function of temperature.At 120 GPa,the free energy value of o-Fe7C3is 14 meV/atom at 0 K and 72 meV/atom at 4000 K, less than the counterparts of h-Fe7C3.At 360 GPa, the free energy of o-Fe7C3is 17 meV/atom more than that of h-Fe7C3at 0 K, but is 29 meV/atom less than that of h-Fe7C3at 4000 K.Thus, the free energy of o-Fe7C3is less than that of h-Fe7C3at high temperatures.In this study,the energy value of self-consistent optimization convergence is set to 10?5eV, which is much smaller than the difference between the free energy values.Therefore, o-Fe7C3appears to be more stable than h-Fe7C3under inner core condition.Hereafter,we investigate the thermal properties of o-Fe7C3(instead of h-Fe7C3[26,27]) at high temperature and high pressure.

    Fig.1.Free energy values of o-Fe7C3 and h-Fe7C3 under simultaneous action of high pressure and high temperature,showing temperature-dependent(a)Gibbs free energy at 120 GPa,(b)Gibbs free energy at 360 GPa,and(c)at 120 GPa temperature-dependent difference in Gibbs free energy between o-Fe7C3 and h-Fe7C3; (d)at 360 GPa temperature-dependent difference in Gibbs free energy between o-Fe7C3 and h-Fe7C3.

    3.2.Thermal EOS of o-Fe7C3

    The Fe7C3is magnetic at low pressures and undergoes magnetic collapse at high pressures.[5]Consequently,we calculate the EOSs for o-Fe7C3with and without spin at 0 K.Figure 2 shows a comparison between the results of this study and previous results from the static compression method at room temperature[5,10,28,29]andab initiosimulations,[13,26,30]whose detailed comparisons can be found in supplementary Table 3.Below 50 GPa, the EOSs with spin consideration are consistent with static measurements at 300 K[5,10,28,29]andab initiosimulations,[13,30]all of which are higher than the simulated results without spin consideration.[13,26]With the increase of pressure, the EOSs with spin restrictions began to converge with the EOSs without spin restrictions due to demagnetization at~70 GPa–100 GPa, implying the presence of a spin phase transition.The transition pressure is consistent with that observed by using the static compression method.[5]The discontinuity in the EOS of Fe7C3near 53 GPa (or 70 GPa) is presumably attributed to a transition from a high-to low-spin state.[5,29]At 300 K, the extrapolations of the experimental measurements are different from each other considerably at inner core pressures(Fig.2),and the results obtained herein accord with the Prescheret al.’s results.[5]The differences from previous experimental results may be related to their samples, some of which were based on powder x-ray diffraction data,[29]whereas the structure reported by Prescheret al.was based on single-crystal x-ray diffraction data.[5]The pressure and density data of the non-magnetic phase fit the third-order Birch–Murnaghan(BM)equation[31]

    wherePis the pressure,ρ0is the initial density,B0is the bulk modulus, andis the pressure derivative.Thus, the bulk modulusB0=322±1 GPa and its pressure derivative=4.38±0.02 can be obtained.

    Figure 3 shows the linear compression ratio of Fe7C3at high pressure.The largest compression ratio is along theaaxis,while the smallest one is along thecaxis,indicating that Fe7C3turns anisotropic under compression.The compressibility along theaaxis is greater than that along thebaxis and thecaxis because the Fe–Fe bond is along theaaxis,which is stronger and tighter than the Fe–C bond along thebaxis andcaxis.Although the compression ratios are different from those observed in static experiments, their variation trends along all three axes are consistent with previous experimental results.[28]

    To ascertain the carbon content of the inner core,we calculate the thermal EOS of o-Fe7C3at high temperature and pressure by using QHA[22](Fig.4).

    Fig.3.Compression ratios varying with pressure in different directions.

    Fig.4.Equation of states of Fe7C3 at specific temperatures from 0 GPa to 360 GPa(dashed lines indicate Fe isotherms).

    In general,the results obtained herein are consistent with the experimental measurements of Lai.[28]In addition,we calculate the EOSs for h-Fe7C3in this study,but they almost coincide with those obtained for o-Fe7C3at high temperature and pressure.The data for pure Fe are also shown here for comparison.[32]At inner core pressures,the density of Fe7C3is 9.47%lower than that of pure Fe at a temperature of 5000 K.

    3.3.Sound velocities

    To ascertain the carbon content of the inner core,it is necessary to determine the elastic properties of Fe7C3at high temperatures and pressures.However, there are clear differences among the experimental measurements of sound velocities in Fe7C3at high pressures (Fig.5), even at room temperature.In this study,the sound velocity of Fe7C3at 0 K is calculated by using the method described in Section 2.Because we are more interested in the sound velocity of Fe7C3at core pressure, we only calculate the sound velocity for Fe7C3without spin restriction(Fig.5).

    Fig.5.Sound velocities of Fe7C3 obtained in this and previous studies,where Xs are cited from the preliminary reference Earth model[PREM],[33] showing density-dependent(a)compressional waves,(b)Debye velocities,and(c)shear wave velocities.

    The calculatedVPof Fe7C3increases linearly with density, complying with Birch’s law:VP=?1.72(±0.12)+1.15(±0.01)ρ.At 0 K, the calculatedVSis also linearly related to the density,VS=0.66(±0.05)+0.34(±0.00)ρ.The Debye velocity(VD)of Fe7C3can be calculated from the following formula:

    The calculatedVP,VD, andVSat 0 K are broadly consistent with the experimental results at 300 K obtained by Prescheret al.,[5]and higher than those obtained by Chen and Zhang,[4]particularly at high pressures.Figure 5 also shows the ab initio values.[27,30]The results obtained herein are similar to those obtained by Liet al.[27]and Daset al.,[34]but they are all less than those calculated by Ghoshet al.[30]Ghoshetal.also suggested that this deviation might be attributed to the nonlinear mixing of phonon modes.[30]

    4.Geophysical implications for inner core composition

    According to the results of this study, the effect of temperature on the sound velocity of Fe7C3is unknown.Our previous experimental results for Fe and Fe–Si alloys suggest that Birch’s law holds true at high temperatures and pressures.[35]Assuming that theVPof Fe7C3follows Birch’s law at high temperatures,we obtain theVPof Fe7C3after determining its density at inner core.The bulk wave velocity(VB)of Fe7C3at high temperatures is calculated according to its EOS given as follows:

    whereKSandKTare the bulk modulus along the isentrope and the isotherm,respectively;ρis the density;αis the thermal expansion coefficient; andγis the Grüneisen parameter.The value ofα,γ, andKTeachare obtained from the phonon spectrum, andVSat high temperatures is obtained from the following formula:

    and Poisson’s ratio is calculated from the following equation:

    The temperature profile of the inner core is calculated by usingT=TICB(ρ/ρICB)γ.The preferred value ofTICBandγare 5400 K and 1.5.[1]However,in this study,we only calculate the data forTICB=5000 K,as QHA often introduces large errors when the temperature is sufficiently high or close to the melting temperature of Fe7C3.[36]The reliability of QHA at high temperatures has been discussed previously.[37]Figure 6 shows the calculated value ofρ,VP,VS, andσof Fe7C3and the comparison with the values of preliminary reference Earth model (PREM).[33]Except for the density, the simulatedVP,VS,andσof Fe7C3obtained in this study deviated from those obtained by Daset al.[34]

    The density of Fe and Fe7C3are about 5.87%higher and 6.34%lower than the PREM values respectively in the whole range of Earth’s inner core pressure.The Fe7C3content of 44.71 wt% can explain the density deficit of the inner core(Fig.6(a))if carbon is the only light element in the core and is present as Fe7C3.Under the inner core condition, the value ofVSandσof Fe7C3are very close to the PREM values compared with those of Fe (Fig.6(c)).TheVPof Fe7C3is slightly higher than that of Fe,and both are much larger than the PREM values.Thus, carbon might be ruled out as the major light element in the Earth’s inner core,which is consistent with findings of Fe–C phase,[38]Fe–C solid solution,[39]and geochemical[40]investigations.By analyzing the Fe–C phase diagram at high temperatures and pressures, Feiet al.suggested that the carbon content of the inner core is about 2.24%.[38]More recently, Huanget al.experimentally determined Fe–C solid solution densities up to core pressures,and concluded that interstitial carbon can lower the density of iron and be present in the inner core.[39]Using geochemical methods,McDonough concluded that the carbon content of the inner core is less than 0.2%.[40]Therefore, if carbon is present in the inner core, other light elements must also be present,which has inspired many studies of the Fe–C–Si ternary system.Pamatoet al.proposed an Fe–C–Si model for the inner core composition based on the EOS of Fe–C–Si at high pressures and temperatures.[41]Daset al.found that doping a small quantity of Si impurities at carbon sites in Fe7C3carbide phases can reduce theVPandVS,which are closer to the PREM data.[34]According to theVPandVSvelocities of Fe and Fe-8.6Si,and combining theVPandVSresults of Fe7C3obtained by Chen and Zhang,[4]Huanget al.[35]suggested that the density and sound velocity for an Fe–C–Si inner core both match the PREM values.However, theVPandVSdata obtained by Chenet al.are less than our calculated results and theab initiovalues obtained by Liet al.,[27]as well as the experimental data obtained by Prescheret al.[5]Therefore, further experiments are required to investigate theVPandVSof Fe7C3,particularly at high temperatures.

    Fig.6.(a)Pressure-dependent calculated densities of Fe(black line)[35] and Fe7C3 (red line), assuming a temperature of 5000 K,compared with the PREM[33] density profile.Fe7C3 (pink line) and Fe7(C,Si)3 (blue dotted line) are cited from Ref.[34].Fe–1Si–5C (green dotted line)is cited from Ref.[35].(b)Pressure-dependent compressional wave velocities and their comparisons.(c)Pressure-dependent shear wave celocities and their comparisons.(d)Pressure-dependent Poisson’s ratios and their comparisons.

    5.Conclusions

    In the present study, the thermal equation of state(EOS)and sound velocities of Fe7C3at high pressures are calculated by using first-principles methods.The simulated results are generally in agreement with the previous experimental results.Although the derivedVSandσof Fe7C3can reproduce the inner core values,the calculatedVPof Fe7C3is higher than that of Fe and the inner core.Therefore,it is less likely for carbon to be present only as Fe7C3in the Earth’s core.

    Acknowledgement

    Project supported by the National Natural Science Foundation of China (Grant Nos.41904085, 41874103, and 42274124).

    猜你喜歡
    莉莉海軍
    曉褐蜻
    綠色天府(2022年6期)2022-07-14 11:59:42
    誰在悄悄幫助莉莉呢?
    我的海軍之夢
    軍事文摘(2020年22期)2021-01-04 02:17:24
    With you at that moment
    相信愛
    不倒自行車
    Look from the Anglo—American jury system of jury system in our country
    封面人物·楊海軍
    新聞愛好者(2016年3期)2016-12-01 06:04:24
    海軍協(xié)議1
    New Approach to Calculate the Unavailability of Identical Spares in Cold Spare Configuration
    亚洲五月婷婷丁香| 日韩免费av在线播放| 久热这里只有精品99| 天天添夜夜摸| 免费在线观看黄色视频的| av网站在线播放免费| 丝袜美足系列| 神马国产精品三级电影在线观看 | 中文字幕色久视频| 久久久久久久久中文| 亚洲一区二区三区色噜噜 | 久久精品国产99精品国产亚洲性色 | 欧美人与性动交α欧美精品济南到| 成人亚洲精品av一区二区 | 99riav亚洲国产免费| 美女扒开内裤让男人捅视频| 亚洲av五月六月丁香网| 国产蜜桃级精品一区二区三区| 丝袜美腿诱惑在线| 中文字幕色久视频| 成年版毛片免费区| 免费av毛片视频| 天堂中文最新版在线下载| 国产成人一区二区三区免费视频网站| 啦啦啦在线免费观看视频4| 久久精品91蜜桃| 精品人妻在线不人妻| 亚洲国产精品一区二区三区在线| 婷婷丁香在线五月| 国产乱人伦免费视频| 在线观看一区二区三区| 电影成人av| 男人操女人黄网站| 欧洲精品卡2卡3卡4卡5卡区| 亚洲午夜理论影院| 久久国产精品影院| 亚洲av美国av| 视频区图区小说| 99精品在免费线老司机午夜| 12—13女人毛片做爰片一| 岛国视频午夜一区免费看| 夜夜爽天天搞| 最近最新中文字幕大全电影3 | 制服人妻中文乱码| 亚洲一区二区三区色噜噜 | 免费av中文字幕在线| 激情在线观看视频在线高清| 他把我摸到了高潮在线观看| 国产成人影院久久av| 国产精品 欧美亚洲| 午夜激情av网站| 男女下面进入的视频免费午夜 | 久久国产精品男人的天堂亚洲| 精品国产乱码久久久久久男人| 极品人妻少妇av视频| 看黄色毛片网站| 18禁观看日本| 91av网站免费观看| 亚洲男人的天堂狠狠| 岛国在线观看网站| 欧美久久黑人一区二区| 亚洲精品中文字幕在线视频| 午夜福利免费观看在线| 亚洲一区高清亚洲精品| netflix在线观看网站| www.自偷自拍.com| 男女之事视频高清在线观看| 欧洲精品卡2卡3卡4卡5卡区| 亚洲中文字幕日韩| 亚洲精品粉嫩美女一区| 国产国语露脸激情在线看| 国产av一区在线观看免费| av免费在线观看网站| 国产又爽黄色视频| 亚洲人成电影观看| 日日夜夜操网爽| 欧美日韩av久久| 久久精品aⅴ一区二区三区四区| 日韩av在线大香蕉| 最新美女视频免费是黄的| 中文亚洲av片在线观看爽| 亚洲成人国产一区在线观看| 免费在线观看黄色视频的| 美女福利国产在线| 男女午夜视频在线观看| 嫩草影院精品99| 日韩一卡2卡3卡4卡2021年| 欧美日韩国产mv在线观看视频| 午夜影院日韩av| 91av网站免费观看| 久久午夜综合久久蜜桃| 亚洲中文字幕日韩| 久久久精品欧美日韩精品| 69精品国产乱码久久久| 欧美人与性动交α欧美软件| 男人舔女人下体高潮全视频| 69av精品久久久久久| 国产单亲对白刺激| 中文亚洲av片在线观看爽| 亚洲成av片中文字幕在线观看| 久久午夜亚洲精品久久| 9191精品国产免费久久| 国产真人三级小视频在线观看| 免费搜索国产男女视频| 午夜免费鲁丝| 久久人人精品亚洲av| 国产伦人伦偷精品视频| 最新美女视频免费是黄的| 久久婷婷成人综合色麻豆| 国产黄a三级三级三级人| 国产国语露脸激情在线看| 91九色精品人成在线观看| 成人三级黄色视频| 这个男人来自地球电影免费观看| 日日夜夜操网爽| 精品人妻1区二区| 午夜亚洲福利在线播放| 丰满迷人的少妇在线观看| 999精品在线视频| 成年人免费黄色播放视频| 精品久久久久久电影网| 岛国在线观看网站| 老汉色av国产亚洲站长工具| 人人妻人人澡人人看| 美女福利国产在线| 欧美黑人精品巨大| 麻豆国产av国片精品| 啦啦啦在线免费观看视频4| 国产男靠女视频免费网站| 天天添夜夜摸| 午夜成年电影在线免费观看| 丝袜在线中文字幕| 亚洲av日韩精品久久久久久密| 久久这里只有精品19| 午夜日韩欧美国产| 国产精品香港三级国产av潘金莲| 超色免费av| 久久热在线av| 国产三级在线视频| 亚洲五月天丁香| 黄片播放在线免费| 久久国产乱子伦精品免费另类| 在线看a的网站| 成人18禁在线播放| 纯流量卡能插随身wifi吗| 精品国产乱码久久久久久男人| 午夜久久久在线观看| 亚洲欧美日韩另类电影网站| 国产一区二区激情短视频| 少妇的丰满在线观看| 老司机在亚洲福利影院| 深夜精品福利| 成人国语在线视频| 亚洲精品在线美女| 性色av乱码一区二区三区2| 高清在线国产一区| 欧美日韩亚洲高清精品| 高清黄色对白视频在线免费看| 国内久久婷婷六月综合欲色啪| 老汉色∧v一级毛片| 俄罗斯特黄特色一大片| 精品久久久久久久久久免费视频 | 国产精品偷伦视频观看了| 国产99久久九九免费精品| 久久人人精品亚洲av| 一边摸一边抽搐一进一小说| 真人做人爱边吃奶动态| 精品人妻在线不人妻| 久久久国产成人免费| 咕卡用的链子| 美女高潮喷水抽搐中文字幕| 可以免费在线观看a视频的电影网站| 亚洲五月色婷婷综合| 亚洲熟妇熟女久久| 午夜久久久在线观看| 免费在线观看影片大全网站| 最近最新中文字幕大全免费视频| 亚洲色图综合在线观看| 男女下面进入的视频免费午夜 | 五月开心婷婷网| 老司机午夜福利在线观看视频| 少妇 在线观看| 成年人免费黄色播放视频| www国产在线视频色| 黑人巨大精品欧美一区二区mp4| 日韩三级视频一区二区三区| 欧美成人免费av一区二区三区| 久久国产精品影院| 久久人妻熟女aⅴ| 久久国产精品人妻蜜桃| 在线观看一区二区三区激情| 国产精品成人在线| www日本在线高清视频| 久久中文看片网| 99久久精品国产亚洲精品| 一二三四在线观看免费中文在| 长腿黑丝高跟| 国产成人啪精品午夜网站| 亚洲国产欧美日韩在线播放| 精品国产乱码久久久久久男人| 久久国产精品人妻蜜桃| 午夜a级毛片| 国产精品美女特级片免费视频播放器 | 久久香蕉精品热| 成在线人永久免费视频| 精品久久久久久久久久免费视频 | 女警被强在线播放| 欧美日韩精品网址| 国产色视频综合| 亚洲欧美一区二区三区久久| 午夜福利,免费看| 亚洲成a人片在线一区二区| 色老头精品视频在线观看| 国产一卡二卡三卡精品| 亚洲成a人片在线一区二区| 亚洲成人免费电影在线观看| 亚洲国产精品一区二区三区在线| 亚洲国产欧美网| 午夜福利一区二区在线看| 日韩有码中文字幕| 亚洲成人免费电影在线观看| 精品福利永久在线观看| 可以在线观看毛片的网站| 国产真人三级小视频在线观看| 久久久国产成人免费| 欧美不卡视频在线免费观看 | 一边摸一边抽搐一进一出视频| 国产欧美日韩精品亚洲av| 亚洲一区高清亚洲精品| 欧美+亚洲+日韩+国产| cao死你这个sao货| 一本大道久久a久久精品| 亚洲熟妇中文字幕五十中出 | 精品乱码久久久久久99久播| 美女扒开内裤让男人捅视频| 欧美av亚洲av综合av国产av| 在线观看日韩欧美| 成人特级黄色片久久久久久久| 99国产综合亚洲精品| 成人18禁高潮啪啪吃奶动态图| 最好的美女福利视频网| 一级,二级,三级黄色视频| 精品高清国产在线一区| 99国产精品免费福利视频| 婷婷丁香在线五月| 色在线成人网| 老司机午夜十八禁免费视频| 亚洲精品成人av观看孕妇| 日韩视频一区二区在线观看| 国产激情欧美一区二区| a级片在线免费高清观看视频| 亚洲一区高清亚洲精品| 好男人电影高清在线观看| 一级,二级,三级黄色视频| 国产成+人综合+亚洲专区| 欧美黑人欧美精品刺激| 欧美大码av| 美女大奶头视频| 久99久视频精品免费| 亚洲五月色婷婷综合| 男女之事视频高清在线观看| 精品久久久精品久久久| 久久久国产精品麻豆| 久久国产精品人妻蜜桃| 免费av毛片视频| 狂野欧美激情性xxxx| 男人舔女人的私密视频| 久久午夜亚洲精品久久| 天堂俺去俺来也www色官网| 天堂中文最新版在线下载| 国产成人免费无遮挡视频| 亚洲成人国产一区在线观看| 高清在线国产一区| 日韩欧美一区视频在线观看| 国产成人欧美在线观看| 韩国精品一区二区三区| 欧美日韩福利视频一区二区| 亚洲一区二区三区色噜噜 | 欧美人与性动交α欧美软件| 亚洲性夜色夜夜综合| 脱女人内裤的视频| 中文字幕人妻丝袜一区二区| 正在播放国产对白刺激| 亚洲片人在线观看| 欧美丝袜亚洲另类 | 国产精品日韩av在线免费观看 | 不卡一级毛片| 18禁裸乳无遮挡免费网站照片 | 欧美另类亚洲清纯唯美| 国产xxxxx性猛交| 国产麻豆69| 久久精品91蜜桃| 一级片免费观看大全| 一本综合久久免费| 日韩视频一区二区在线观看| 亚洲av片天天在线观看| 午夜福利免费观看在线| 男男h啪啪无遮挡| 岛国视频午夜一区免费看| 日本三级黄在线观看| 在线视频色国产色| 夜夜夜夜夜久久久久| 老司机靠b影院| 1024视频免费在线观看| 欧美不卡视频在线免费观看 | 午夜两性在线视频| 99riav亚洲国产免费| 成年版毛片免费区| 国产精品久久久久久人妻精品电影| 精品国产国语对白av| 啦啦啦 在线观看视频| 中文字幕最新亚洲高清| 久久精品国产清高在天天线| 国产精品永久免费网站| 国产av一区二区精品久久| 久久草成人影院| 国产深夜福利视频在线观看| 人人妻人人澡人人看| 悠悠久久av| 18禁裸乳无遮挡免费网站照片 | 成人18禁在线播放| 日本 av在线| x7x7x7水蜜桃| 新久久久久国产一级毛片| 在线观看日韩欧美| 男人舔女人的私密视频| 久久久精品国产亚洲av高清涩受| 韩国精品一区二区三区| 很黄的视频免费| 国产精品综合久久久久久久免费 | 性少妇av在线| 亚洲精品在线观看二区| 很黄的视频免费| 自线自在国产av| 国产欧美日韩一区二区精品| 国产成人欧美| 久久久久九九精品影院| 国产亚洲精品久久久久久毛片| 交换朋友夫妻互换小说| 9热在线视频观看99| 香蕉丝袜av| 最近最新中文字幕大全电影3 | 亚洲自偷自拍图片 自拍| av在线天堂中文字幕 | 精品乱码久久久久久99久播| 精品久久久久久久久久免费视频 | 亚洲男人的天堂狠狠| 日韩 欧美 亚洲 中文字幕| 女人被狂操c到高潮| 日本五十路高清| 制服诱惑二区| 亚洲成人久久性| 欧美色视频一区免费| 亚洲一区二区三区色噜噜 | 国产精品一区二区精品视频观看| 啪啪无遮挡十八禁网站| 欧美国产精品va在线观看不卡| 露出奶头的视频| 免费在线观看日本一区| 免费高清视频大片| 人妻丰满熟妇av一区二区三区| 99热只有精品国产| 精品国产乱子伦一区二区三区| 亚洲成av片中文字幕在线观看| 国产精品一区二区免费欧美| 每晚都被弄得嗷嗷叫到高潮| 久久久国产成人免费| 麻豆久久精品国产亚洲av | 亚洲欧美日韩无卡精品| 久久久久久亚洲精品国产蜜桃av| 亚洲国产毛片av蜜桃av| 91麻豆av在线| 久久精品国产亚洲av香蕉五月| 日本三级黄在线观看| 两个人免费观看高清视频| 日韩一卡2卡3卡4卡2021年| 精品人妻1区二区| 操美女的视频在线观看| 新久久久久国产一级毛片| 免费不卡黄色视频| 久久久国产欧美日韩av| 国产亚洲欧美精品永久| 国产精品香港三级国产av潘金莲| 久久久久精品国产欧美久久久| 亚洲视频免费观看视频| 精品国产乱子伦一区二区三区| 亚洲精品美女久久av网站| 国产激情欧美一区二区| 国产一区二区三区在线臀色熟女 | 交换朋友夫妻互换小说| 在线观看免费视频日本深夜| 国产精华一区二区三区| 亚洲精品国产区一区二| 人人妻人人澡人人看| 99国产极品粉嫩在线观看| 国产亚洲欧美98| 久久人人精品亚洲av| 中出人妻视频一区二区| 黄色女人牲交| 两个人看的免费小视频| 久久久水蜜桃国产精品网| 窝窝影院91人妻| 19禁男女啪啪无遮挡网站| 一区二区三区国产精品乱码| 俄罗斯特黄特色一大片| 亚洲午夜精品一区,二区,三区| 亚洲av第一区精品v没综合| 免费日韩欧美在线观看| 视频在线观看一区二区三区| 国产无遮挡羞羞视频在线观看| 成人黄色视频免费在线看| 午夜精品国产一区二区电影| 国产97色在线日韩免费| 一a级毛片在线观看| 国产精品 欧美亚洲| 国产精品免费一区二区三区在线| 久久久久久久久免费视频了| 国产一区二区在线av高清观看| 一级毛片精品| 91av网站免费观看| 欧美精品亚洲一区二区| 黄片大片在线免费观看| 99久久久亚洲精品蜜臀av| 三级毛片av免费| 午夜影院日韩av| av网站在线播放免费| 丰满的人妻完整版| 亚洲 欧美一区二区三区| 一个人免费在线观看的高清视频| 亚洲专区字幕在线| 日韩av在线大香蕉| 热re99久久精品国产66热6| 免费在线观看黄色视频的| 国产男靠女视频免费网站| 成年女人毛片免费观看观看9| 亚洲片人在线观看| 亚洲欧美日韩高清在线视频| 成在线人永久免费视频| 亚洲色图av天堂| a在线观看视频网站| 每晚都被弄得嗷嗷叫到高潮| 国产97色在线日韩免费| 精品一区二区三卡| 最新美女视频免费是黄的| 日韩三级视频一区二区三区| 国产成+人综合+亚洲专区| 免费搜索国产男女视频| 国产一区二区三区视频了| 黄色片一级片一级黄色片| 亚洲精品一区av在线观看| 久久精品国产99精品国产亚洲性色 | 国产激情欧美一区二区| 午夜福利免费观看在线| 激情在线观看视频在线高清| 久久午夜亚洲精品久久| 最好的美女福利视频网| 国内久久婷婷六月综合欲色啪| 男女床上黄色一级片免费看| 91在线观看av| 久久久久久久久久久久大奶| 丝袜在线中文字幕| 欧美中文日本在线观看视频| 国产欧美日韩一区二区三| 亚洲 欧美一区二区三区| 国产av精品麻豆| 日韩精品中文字幕看吧| 一区二区三区精品91| 自线自在国产av| 欧美一级毛片孕妇| 深夜精品福利| 香蕉久久夜色| 日韩大码丰满熟妇| 欧美+亚洲+日韩+国产| 精品欧美一区二区三区在线| 99热只有精品国产| 精品久久久精品久久久| 视频在线观看一区二区三区| 亚洲全国av大片| 性欧美人与动物交配| av视频免费观看在线观看| 日日夜夜操网爽| 久久人妻熟女aⅴ| 中文亚洲av片在线观看爽| 人人妻人人澡人人看| 国产一区在线观看成人免费| 淫妇啪啪啪对白视频| 午夜福利免费观看在线| 777久久人妻少妇嫩草av网站| 成人黄色视频免费在线看| 亚洲美女黄片视频| 国产片内射在线| 女性生殖器流出的白浆| 老熟妇乱子伦视频在线观看| 久久香蕉激情| 国产真人三级小视频在线观看| 亚洲情色 制服丝袜| 黄色视频,在线免费观看| 亚洲专区国产一区二区| 热99re8久久精品国产| 一级a爱片免费观看的视频| 久久热在线av| 香蕉丝袜av| 国产亚洲精品久久久久久毛片| 久久久精品国产亚洲av高清涩受| 精品人妻在线不人妻| 亚洲精品av麻豆狂野| 午夜亚洲福利在线播放| 久久久久久免费高清国产稀缺| 极品人妻少妇av视频| 午夜福利一区二区在线看| 久久国产亚洲av麻豆专区| 99国产精品一区二区蜜桃av| 亚洲一区二区三区色噜噜 | 嫩草影视91久久| 国产精品久久电影中文字幕| av有码第一页| 9热在线视频观看99| 中文字幕人妻熟女乱码| 国产真人三级小视频在线观看| 夜夜躁狠狠躁天天躁| 日韩欧美免费精品| 999久久久精品免费观看国产| 两人在一起打扑克的视频| 欧美午夜高清在线| 色播在线永久视频| 麻豆av在线久日| 久久欧美精品欧美久久欧美| 精品久久久久久成人av| 日韩三级视频一区二区三区| 手机成人av网站| 在线观看日韩欧美| 亚洲精品久久成人aⅴ小说| 久久 成人 亚洲| 国产精品久久久久久人妻精品电影| 久久久久久人人人人人| www.精华液| 色综合站精品国产| 一级黄色大片毛片| 亚洲成a人片在线一区二区| 国产一卡二卡三卡精品| 成人永久免费在线观看视频| 超色免费av| 波多野结衣av一区二区av| 最近最新中文字幕大全电影3 | 亚洲第一av免费看| 两个人看的免费小视频| 天堂中文最新版在线下载| 色哟哟哟哟哟哟| 欧美中文日本在线观看视频| 精品一品国产午夜福利视频| 国产aⅴ精品一区二区三区波| www.999成人在线观看| 国产免费男女视频| 久久精品91无色码中文字幕| 日本黄色视频三级网站网址| 99热国产这里只有精品6| 在线观看一区二区三区激情| 脱女人内裤的视频| 夜夜爽天天搞| av天堂在线播放| 丰满迷人的少妇在线观看| 午夜亚洲福利在线播放| 12—13女人毛片做爰片一| 精品久久久久久成人av| 露出奶头的视频| 亚洲熟女毛片儿| 十分钟在线观看高清视频www| 亚洲av美国av| 亚洲国产精品合色在线| 黄色片一级片一级黄色片| 欧美午夜高清在线| 久久精品亚洲精品国产色婷小说| 男人操女人黄网站| 国产主播在线观看一区二区| 精品高清国产在线一区| 亚洲人成电影观看| 一级a爱片免费观看的视频| 丁香六月欧美| 国产黄色免费在线视频| 亚洲美女黄片视频| 欧美最黄视频在线播放免费 | 久久精品亚洲精品国产色婷小说| 男人操女人黄网站| 日韩欧美一区视频在线观看| 咕卡用的链子| 欧美色视频一区免费| 男人操女人黄网站| 在线观看免费视频网站a站| 精品卡一卡二卡四卡免费| 成人手机av| 啦啦啦免费观看视频1| 成人精品一区二区免费| 黄片大片在线免费观看| 99re在线观看精品视频| 亚洲av成人一区二区三| 一区二区三区精品91| 国产又色又爽无遮挡免费看| 精品卡一卡二卡四卡免费| 成年人黄色毛片网站| 精品国产一区二区三区四区第35| 国产1区2区3区精品| 午夜免费激情av| 国产精品久久久av美女十八| 久久国产亚洲av麻豆专区| 嫁个100分男人电影在线观看| 亚洲国产精品999在线| 真人做人爱边吃奶动态| 午夜免费激情av| 久久狼人影院| 咕卡用的链子| 别揉我奶头~嗯~啊~动态视频| 精品国产美女av久久久久小说| 日韩国内少妇激情av| 日本wwww免费看| 亚洲av成人av|