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

雙磁性中心內嵌富勒烯Y2C2@C82-C2(1)中的超快自旋動力學行為*

2019-03-11 08:54:44黃瑞李春金蔚GeorgiosLefkidisWolfgangbner西北工業大學力學與土木建筑學院西安71007陜西師范大學物理學與信息技術學院西安710119德國凱澤斯勞滕工業大學物理系OPTIMAS研究中心凱澤斯勞滕67653德國2018年10月23日收到2018年12月7日收到修改稿
物理學報 2019年2期
關鍵詞:理論體系

黃瑞 李春? 金蔚 Georgios Lefkidis Wolfgang Hübner (西北工業大學力學與土木建筑學院, 西安 71007 (陜西師范大學物理學與信息技術學院, 西安 710119) (德國凱澤斯勞滕工業大學物理系, OPTIMAS研究中心, 凱澤斯勞滕 67653, 德國)(2018 年10 月23日收到; 2018 年12 月7日收到修改稿)

自旋翻轉和自旋轉移是實現基于內嵌富勒體系自旋邏輯功能器件設計的先決條件. 本文以雙磁性中心內嵌富勒烯Y2C2@C82-C2(1)體系為例, 采用第一性原理計算方法, 結合Λ進程理論模型和自編的遺傳算法程序, 在該內嵌富勒烯體系中分別實現了亞皮秒時間尺度內的自旋翻轉和自旋轉移過程. 計算結果表明, 優化后的內嵌Y2C2團簇結構和實驗得到的各項數據基本吻合, 并且會對外部的C82-C2(1)籠結構產生一定的排斥力, 但由于富勒烯籠狀結構具有很強的穩定性, 所以整個體系仍然保持碳籠結構的完整性. 通過對自旋密度分布與激光脈沖作用下自旋期望值演化的具體分析, 經由Λ進程的自旋翻轉是基于兩個Y元素的整體自旋翻轉; 自旋轉移則源自兩個磁性中心以及碳籠之間在激光脈沖作用下的自旋密度重新分布. 本文結果揭示了Y2C2@C82-C2(1)體系中的超快自旋動力學機理, 可望為基于實際內嵌富勒烯分子的自旋邏輯功能器件設計提供理論依據.

1 引 言

內嵌富勒烯由于其獨特的空心籠狀結構可以內嵌入不同的原子、分子、離子或團簇, 使其不但具有富勒烯碳籠的物理化學性質, 而且還兼具了內嵌團簇的各種優良特性, 在磁學[1-3]、力學[4,5]、化學[6,7]、半導體材料[8,9]和生物醫學[10,11]等領域擁有巨大的潛在應用價值, 被認為是未來納尺度功能器件中重要的備選材料之一. 1991年, 通過激光蒸發氧化鑭/石墨復合棒, 第一個宏觀量的內嵌富勒烯La@C82在實驗中被合成出來[12], 進而得到了內嵌富勒烯配合物的詳盡波譜. 隨著科研人員對富勒烯配合物研究的不斷深入, 越來越多新內嵌富勒烯結構被發現和實驗分離出來[13-17], 其內嵌化合物也從最開始的單個原子發展到分子、離子甚至分子團簇. 內嵌富勒烯也由于其內嵌團簇的特殊電磁特性, 在納米級邏輯功能器件中具有廣泛的潛在應用前景. 人們發現內嵌富勒烯可以作為量子信息的載體, 通過對其電荷和自旋的操控, 設計出可擴展量子計算機體系[18-20], 從而有望大幅度提高計算設備的運算速度和存儲密度.

目前, 基于電荷自由度的傳統電子學仍在遵循著摩爾定律發展. 然而隨著集成電路上線路密度的增加, 其差錯率和復雜性也將呈指數增長. 由于人們對于高性能計算條件的需求不斷提高, 必須考慮采用新技術和新材料來滿足不斷增長的需求, 因此人們將目光聚集到通過操控自旋自由度來提高計算器件的運算速度和讀寫速度等關鍵性能上.1996年, 法國科學家Beaurepaire等[21]首次發現了鐵磁金屬中由激光誘導的超快磁化動力學行為, 即通過實驗實現了激光誘導的Ni薄膜中的超快退磁過程, 這一現象的發現馬上引起了國際上研究人員對于各種光驅動的磁化操控機理的研究熱情.2005年, Koopmans團隊[22]提出了聲子調制Elliott-Yafet機制, 利用磁光抽運探測技術研究了鐵磁鎳中電荷和自旋的弛豫過程. 2009年, Bigot研究組[23]通過研究激光脈沖與鐵磁金屬中的電子自旋耦合作用, 提出了相對論量子電動力學理論. 2010年, Battiato等[24]提出了由飛秒激光誘導退磁導致的自旋極化激發電子的超擴散機制. St?hr和他的同事[25]通過大量的理論研究, 提出了Λ進程理論模型, 用來描述激光與磁性體系的相互作用, 以及激光誘導的超快光學躍遷過程. 通過對這一理論進行修正與完善, Li等[5,26,27]應用第一性原理方法實現了內嵌磁性元素富勒烯中亞皮秒時間尺度的超快自旋翻轉和自旋轉移過程. 然而, 之前的研究內容主要涉及內嵌富勒烯的理論預測體系, 并未將現有理論模型應用于實驗中已合成的實際內嵌富勒烯分子.

本文選取雙磁性中心內嵌富勒烯Y2C2@C82-C2(1)體系[28]為研究對象, 采用基于Λ進程理論模型的第一性原理計算方法和自編程序, 對該分子體系在磁光耦合作用下的超快自旋動力學行為進行了系統研究, 并實現了在同一種體系中的超快自旋翻轉和自旋轉移. 本文結果可望為未來納尺度磁性存儲單元和自旋邏輯器件設計提供理論依據和新的思路.

2 理論方法

本文實現的超快自旋翻轉和自旋轉移是基于Λ進程理論模型所得到的. 該理論模型, 即由相干激光場驅動的發生在磁性材料中經由中間態參與的初末態之間的間接光學躍遷機理, 已經在我們前期研究工作中進行了詳細的描述[29-32]. 其中, 當選取從同一個三重態通過Zeeman作用分裂為具有相反自旋方向的兩個簡并態分別作為初始態和最終態時, 在特定激光脈沖作用后, 若體系從初始態激發至自旋混合的中間態, 再由中間態返回至最終態, 即認為實現了體系自旋方向的翻轉過程, 其翻轉速度可以達到亞皮秒級別. 另外, 當初始態和最終態的自旋密度高度局域化在不同的磁性原子上時, 則Λ進程也可以實現體系的自旋轉移過程. 具體理論計算方法介紹如下.

首先, 在無任何外場作用時, 可以利用量子化學計算軟件Gaussian 09[33]軟件包對體系的非相對論哈密頓量進行求解, 體系的哈密頓量為

(1)式中的Nat和Nel分別表示原子和電子的數目;Ra和 ri分別為原子和電子的坐標矢量;Za為核電荷數. 由于磁性體系較強的電子相互作用, 要得到較為精確的多體電子基態和激發態, 需要采用高精度的post-Hartree-Fock方法, 本文計算時采用了對稱性匹配簇-組態相互作用方法(symmetry adapted cluster-configuration interaction, SACCI)[34]. 在Λ進程理論模型中, 考慮到體系自旋轉換過程中初始態、中間態和最終態的選取, 需要對體系施加一個合適的外加靜磁場和加入自旋軌道耦合效應(spin-orbit coupling, SOC)來實現Zeeman分裂, 進而產生自旋混合態, 相應的哈密頓量為

3 結果與討論

3.1 Y2C2@C82-C2(1)結構分析

目前, 對于內嵌富勒烯中的超快自旋動力學行為研究大部分僅限于理論預測的單磁性中心內嵌富勒烯體系. 理論上來說, 自旋轉移和自旋翻轉是實現內嵌富勒體系自旋邏輯功能器件設計必不可少的先決條件. 單磁性中心內嵌富勒烯體系可以實現其自旋翻轉操控, 由于這種體系只具有一個磁性中心, 所以單磁性中心體系是無法實現自旋轉移的. 而自旋邏輯器件設計中很重要的一個方面是可以將自旋翻轉、自旋轉移以及各種輸入輸出參數相結合, 以實現相對復雜的邏輯運算功能. 本文選取雙磁性中心內嵌富勒烯Y2C2@C82-C2(1)作為研究對象, 基于Λ進程理論模型和第一性原理計算, 并在前期研究工作的基礎上, 集中探討該分子體系在磁光耦合作用下的自旋動力學行為.

在體系結構選擇上, 首先對C82碳籠結構的9種異構體C2(1),Cs(2),C2(3),Cs(4),C2(5),Cs(6),C3v(7),C3v(8),C2v(9)進行結構優化計算并得到相應的穩定結構, 再將蝶形雙金屬碳化物團簇(Y2C2)內嵌到優化好的C82結構中, 然后再次進行結構優化計算得到了7種穩定的Y2C2@C82異構體C2(1),Cs(2),C2(3),Cs(4),Cs(6),C3v(8),C2v(9), 但是在進行自旋密度計算時只有C2(1)構型中自旋密度是高度局域化在兩個Y原子上, 即最適合本文采用的Λ進程理論模型進行分析求解,因此本文選取的研究對象為Y2C2@C82-C2(1)異構體. 在計算過程中, 為了保證計算效率和計算精度,在選取基組時, 對于電子結構相對簡單的C原子,選取STO-3G基組, 對于電子結構相對較為復雜的Y原子, 選取帶有有效核勢(ECP)的LANL2DZ基組[35]. 本文計算中所采用的收斂標準為: 計算得到的每個原子上的力和力的均方根分別小于4.5 ×10-4和3 × 10-4Hartree/Bohr, 其相應的位移和位移的均方根則分別小于 1.8 × 10-3和 1.2 ×10-3Bohr.

圖1 雙磁性中心內嵌富勒烯Y2C2@C82-C2(1)體系 (a) 內嵌富勒烯Y2C2@C82C2(1); (b) 內嵌的雙金屬碳化物團簇(Y2C2)Fig.1. Double-magnetic-center endohedral fullerene Y2C2@C82-C2(1) system: (a) Endohedral fullerene Y2C2@C82C2(1);(b) endohedral bimetallic carbide clusters (Y2C2).

3.2 Y2C2@C82-C2(1)中的超快自旋翻轉過程

內嵌金屬富勒烯作為自旋邏輯功能器件潛在材料的最大優勢是其外部的富勒烯籠狀結構可以為內嵌磁性團簇提供一個較為穩定的化學環境, 從而保證內嵌磁性團簇的自旋狀態不受外部環境因素的干擾. 但是在實際應用中, 內嵌磁性團簇有可能會與外部富勒烯籠結構發生相互作用, 使自旋密度從磁性團簇分散到外部富勒烯籠結構上, 使得內嵌磁性團簇的自旋密度局域化程度降低, 導致對其自旋操控變得更為困難. 因此, 需要首先對內嵌富勒烯Y2C2@C82-C2(1)進行自旋密度分布計算. 本文采用SAC-CI作用方法, 計算得到了該體系中各個原子上的自旋密度分布情況, 如表1所列. 通過表1可以看出, 在Y2C2@C82-C2(1)結構基態三重態中, 整個體系中的自旋密度極少分布在C原子上, 而三重態|1〉,|2〉和|4〉上的自旋密度高度局域化在兩個Y原子上, 這樣在選取合適的初始態、中間態和最終態的前提下, 可以通過Λ進程, 實現在激光誘導作用下該體系中Y原子上的自旋翻轉過程.

依據Λ進程理論, 實現整個體系的自旋翻轉過程最簡單有效的是選擇兩個加入SOC作用之前簡并的、且自旋方向相反的自旋態作為初始態和最終態. 施加SOC和外磁場作用后, 在特定激光脈沖的驅動下, 體系若可以經由中間態完成從初始態到最終態的躍遷, 則實現了自旋翻轉過程. 本文通過對Y2C2@C82-C2(1)體系施加一個磁場強度為2.35 T的外加靜磁場, 使得該體系中所有的三重態發生Zeeman分裂, 第一個三重態分裂得到兩個自旋方向相反的態|2〉和|4〉, 作為自旋翻轉過程的初末態. 在合適激光脈沖作用下, 得到了如圖2所示的自旋翻轉動力學過程, 其中自旋態|2〉(如圖2(d)中黑色虛線)自旋方向向下, 其自旋平均值為-0.962, 自旋態|4〉(如圖2(d)中紅色實線)自旋方向向上, 其自旋平均值為0.963. 整個自旋翻轉過程所涉及到的中間態包括|3〉,|2 9〉,|3 0〉和|3 1〉, 其中主要參與自旋翻轉過程的中間態為能態|2 9〉,|30〉和|3 1〉(如圖2(d)中綠色虛線). 圖2(a)展示了該體系中各能態隨時間變化的占據情況. 可以看出, 在外加激光脈沖作用下時, 初始態|2〉的占據率隨時間逐漸由1降低為0, 而最終態|4〉的占據率則由最初的0上升至最終的0.990, 實現了該體系中Y原子上的自旋狀態由自旋向上轉變為自旋向下的超快自旋翻轉過程, 整個翻轉過程耗時約1000 fs,其保真度高達97.8%. 圖2(b)展示了x,y,z3個坐標軸方向, 體系自旋平均值分量隨時間的變化情況(藍色實線為z軸分量, 紅色實線為y軸分量,黑色實線為x軸分量), 其中分量〈Sz〉主導著整個Λ進程中體系自旋方向的變化狀態. 在此進一步對比體系總自旋平均值的變化情況, 在激光脈沖作用之前, 體系總自旋平均值為0.962, 其方向向量為(-0.001, -0.043, 0.962), 激光作用后,總自旋自旋平均值變為0.963, 其方向向量為(-0.001,-0.040, 0.962), 二者夾角大小為 176.80°, 非常接近180°, 進一步地證實了該體系在激光作用下實現了自旋翻轉過程. 另外, 如圖2(c)所示, 實現超快自旋翻轉過程所采用的激光為線偏振光, 激光能量Elaser為1.59 eV, 半高寬FWHM為300.00 fs,激 光 方 向 在 球 坐 標 系 下 的 入 射 角θ和φ分 別 為142.87°和247.39°, 激光偏振方向與光學平面間的夾角γ為 254.08°.

表1 Y2C2@C82-C2(1)體系中能量最低的五個三重態的自旋密度分布情況(未考慮SOC)Table 1. Spin density distributions of the five lowest triplet states in Y2C2@C82-C2(1) system (without SOC).

表2 自旋翻轉過程涉及能態的局域化位置、能量和自旋期望值(考慮SOC)Table 2. Localization positions, energies and spin expectation values of the involved states in spin-switching scenario (including SOC).

圖2 激光脈沖作用下, Y2C2@C82-C2(1)體系中的超快自旋翻轉過程 (a) 自旋翻轉過程中所涉及的初始態(黑色虛線)、最終態(紅色實線)和中間態(彩色實線)隨時間變化的占據情況; (b) 自旋角動量分量的期望值隨時間的變化情況; (c) 實現自旋翻轉所用的激光脈沖包絡線, 插圖為結構優化后的內嵌Y2C2團簇, 碳籠用虛線表示; (d)SAC-CI計算得到的能級圖(考慮SOC)Fig.2. Ultrafast spin-switching scenario achieved in Y2C2@C82-C2(1) system under laser pulses: (a) Time evolution of the occupation of the initial state (dashed black line), final state (solid red line), and the intermediate states (colorized solid lines)involved in spin-switching scenario; (b) variation of the expectation values of the spin angular momentum components along the x,yand zaxis with time; (c) laser pulse envelope of the spin-switching scenario, inset represents structurally optimized endohedral Y2C2 cluster, the dashed circle represents the carbon cage; (d)energy levels of Y2C2@C82-C2(1) calculated by SAC-CI (including SOC).

3.3 Y2C2@C82-C2(1)中的超快自旋轉移過程

實現自旋轉移是在內嵌富勒烯結構中實現自旋邏輯功能的先決條件之一, 然而最近的理論研究表明, 經由Λ進程自旋動力學理論機理, 自旋轉移相比于自旋翻轉更難實現. 雙磁性中心內嵌富勒烯Y2C2@C82-C2(1)體系擁有兩個磁性中心Y1和Y2原子, 考慮SOC作用后, 得到了該體系中各個原子上自旋密度重新排布情況(如表3所示). 能態|2〉上Y1和Y2原子的自旋密度分別為0.504和0.680 (自旋方向向下), 可以看作對等分布. 單個C原子上自旋密度最大值僅為0.016, 84個C原子上總的自旋密度為1.062, 說明自旋密度高度局域化在Y原子上. 而能態|1 1〉上Y1和Y2原子的自旋密度分別為0.414和0.900 (自旋向下), 單個C原子上自旋密度最大值僅為0.032, 84個C原子上總的自旋密度為0.933, 可以認為自旋密度局域化在Y2原子上. 根據Λ進程理論, 選取能態|2〉作為初始態, 能態|1 1〉作為最終態.

圖3通過Λ進程實現了初始態|2〉到最終態|11〉上的自旋轉移過程. 在外加激光脈沖作用下,能態|2〉的占據率隨時間迅速從1降低為0, 而能態|11〉的占據率則由0上升至0.963, 整個自旋轉移過程耗時約200 fs, 其保真度高達95.1%. 結合表3中自旋密度分布情況, 在整個自旋轉移過程中,Y1原子和碳籠上的自旋密度值均有所減小, 而Y2原子上的自旋密度值明顯增大, 因此, 可以看作是實現了從Y1原子到Y2原子上的自旋轉移過程.值得注意的是, 整個自旋轉移過程涉及到的中間態較 多 (態|1〉,|3〉,|4〉,|12〉,|13〉,|14,|15〉,|16,|29,|30〉和|31〉, 如圖3(d)中綠色虛線所示), 但實際主要參與的中間態僅有態|14〉和|29.另外, 如圖3(c)所示, 實現超快自旋轉移過程所采用的激光為線偏振光, 激光能量Elaser為0.83 eV, 半高寬FWHM為50.75 fs, 激光方向在球坐標系下的入射角θ和φ分別為 253.72°和 305.81°, 激光偏振方向與光學平面間的夾角γ為20.41°.

表3 Y2C2@C82-C2(1)體系中自旋密度分布情況(考慮SOC)Table 3. Spin density distribution in Y2C2@C82-C2(1) system (including SOC).

圖3 激光脈沖作用下, Y2C2@C82-C2(1)體系中的超快自旋轉移過程 (a) 自旋轉移過程中所涉及的初始態(黑色虛線)、最終態(紅色實線)和中間態(彩色實線)隨時間變化的占據情況; (b) 自旋角動量分量的期望值隨時間的變化情況; (c) 實現自旋轉移所用的激光脈沖包絡線, 插圖為結構優化后的內嵌Y2C2團簇, 碳籠用虛線表示; (d)SAC-CI計算得到的能級圖(考慮SOC)Fig.3. Ultrafast spin-transfer scenario achieved in Y2C2@C82-C2(1) system under laser pulses: (a) Time evolution of the occupation of the initial state (dashed black line), final state (solid red line), and the intermediate states (colorized solid lines) involved in spin-transfer scenario; (b) variation of the expectation values of the spin angular momentum components along the x, yand zaxis with time; (c) laser pulse envelope of the spin-transfer scenario, inset represents structurally optimized endohedral Y2C2 cluster, the dashed circle represents the carbon cage; (d) energy levels of Y2C2@C82-C2(1) calculated by SAC-CI (including SOC).

圖4展示了Y2C2@C82-C2(1)體系初始態|2〉和最終態|1 1〉的電子吸收光譜圖, 橫軸表示體系各能態與選取的初始態/最終態之間的能量差值, 縱軸表示振子強度. 電子吸收光譜峰意味著相應態與初始態|2〉/最終態|11〉的躍遷概率不為0, 多峰值情況表示有多個可與初末態發生躍遷的電子態存在. 但是, 考慮到實現Λ進程對兩分支(初態—中間態躍遷與中間態—末態躍遷)對稱性的要求以及激光失諧等因素, 并非所有的峰值所對應的能態都能作為自旋翻轉過程的中間態. 對于自旋翻轉過程, 圖4(a)給出了初始態|2〉的電子吸收光譜(最終態|4〉由于與初始態|2〉從同一個三重態分裂出來,其吸收光譜基本相同, 故沒有顯示). 其中藍色虛線對應的峰值所對應的能態為自旋翻轉過程中涉及到的主要中間態|2 9〉,|30〉和|31〉, 與初始態|2〉之間的能量差值為1.58 eV (其能級位置也可由圖2(d)看出). 該能量差值與實現該自旋翻轉過程所需激光脈沖的能量 (1.59 eV)基本相同, 表明激光失諧較小. 對于自旋轉移過程, 圖4(b)給出了最終態|11〉的電子吸收光譜. 由于選取的初始態|2〉和最終態|1 1〉來自兩個不同的三重態, 其光學性質不盡相同. 經分析發現參與該過程的中間態只要有態|14〉和態|2 9〉, 其位置分別由紅色虛線(圖4(a)能量差0.76 eV處及圖4(b)能量差0.03 eV處)和藍色虛線 (圖 4(a)能量差 1.58 eV處及圖 4(b)能量差0.857 eV處)標注. 其中, 該自旋轉移過程所需外加激光脈沖的能量為0.85 eV, 與E|14〉-|2 9〉能量差基本吻合, 因此易于實現.

通過體系電子吸收光譜可以發現, 自旋翻轉和自旋轉移過程中涉及到的主要中間態與初末態均有躍遷. 因此, 通過分析初始態對應的電子吸收光譜, 可以預判出自旋翻轉和自旋轉移所采用的激光脈沖能量及其實現自旋動力學過程的難易程度.

4 結 論

圖4 Y2C2@C82-C2(1)體系自旋翻轉和自旋轉移過程中初始態 |2〉和最終態 |1 1〉的電子吸收光譜圖Fig.4. Electronic absorption spectra of the initial (state |2 〉) and final (state |1 1〉) state of the spin-switching and spin-transfer scenarios in Y2C2@C82-C2(1).

本文以雙磁性中心內嵌富勒烯Y2C2@C82-C2(1)體系為例, 采用第一性原理計算方法, 結合Λ進程理論模型和自編的遺傳算法程序, 系統地研究了該分子體系在磁光耦合作用下的超快自旋動力學行為. 計算結果表明, 優化后的內嵌Y2C2團簇結構和實驗得到的各項數據基本吻合,并且會對外部的C82-C2(1)籠結構產生一定的排斥力, 但由于富勒烯籠狀結構具有很強的穩定性, 所以整個體系仍然保持碳籠結構的完整性. 在Y2C2@C82-C2(1)體系中, 通過計算發現自旋密度高度局域化在兩個Y原子上, 極少地分布在碳籠上. 經由Λ進程理論模型, 在同一內嵌富勒烯體系中分別實現自旋翻轉和自旋轉移過程. 在自旋翻轉過程中, 本文實現了該體系中Y原子上的自旋狀態由自旋向上轉變為自旋向下的超快自旋翻轉過程, 整個翻轉過程耗時約1000 fs, 其保真度高達97.8%; 同時本文也實現了該體系下的超快自旋轉移過程, 整個自旋轉移過程耗時僅200 fs, 其保真度高達95.1%. 另外, 在進行電子吸收光譜計算過程中發現, 當躍遷所需能量差值與外加激光脈沖能量的失諧程度越小, 實現自旋翻轉和自旋轉移過程的概率越大. 本文的理論研究結果可為基于多磁性中心內嵌富勒烯體系的自旋邏輯功能器件設計奠定理論基礎并提供新的設計思路.

猜你喜歡
理論體系
堅持理論創新
當代陜西(2022年5期)2022-04-19 12:10:18
神秘的混沌理論
理論創新 引領百年
構建體系,舉一反三
相關于撓理論的Baer模
探索自由貿易賬戶體系創新應用
中國外匯(2019年17期)2019-11-16 09:31:14
理論宣講如何答疑解惑
學習月刊(2015年21期)2015-07-11 01:51:44
如何建立長期有效的培訓體系
現代企業(2015年1期)2015-02-28 18:43:18
“曲線運動”知識體系和方法指導
“三位一體”德育教育體系評說
中國火炬(2010年7期)2010-07-25 10:26:09
主站蜘蛛池模板: 亚洲精品卡2卡3卡4卡5卡区| 无码精品国产VA在线观看DVD| 国产乱子伦精品视频| 欧美综合激情| 99热国产这里只有精品无卡顿"| 国产va欧美va在线观看| 无码精品国产dvd在线观看9久| 国产91视频观看| 久一在线视频| 亚洲日本韩在线观看| 99999久久久久久亚洲| 国产免费人成视频网| 免费观看成人久久网免费观看| 国产精品手机在线观看你懂的| 夜夜操天天摸| 午夜视频在线观看免费网站| 日本高清免费不卡视频| 91久久国产热精品免费| 国产精品丝袜视频| 精品伊人久久大香线蕉网站| 亚洲自偷自拍另类小说| 国模视频一区二区| 午夜毛片免费看| 欧美午夜在线视频| yy6080理论大片一级久久| 91色在线观看| 色噜噜在线观看| 亚洲欧洲日韩国产综合在线二区| 亚洲视频在线网| 囯产av无码片毛片一级| 欧美精品一二三区| 免费Aⅴ片在线观看蜜芽Tⅴ| 天天综合网在线| 欧美日韩成人在线观看| 国产导航在线| 亚洲精选高清无码| 五月天福利视频| 欧美视频在线播放观看免费福利资源| 色成人综合| 亚洲精品无码人妻无码| 中文天堂在线视频| 欧美精品高清| 亚洲91在线精品| 在线看片中文字幕| 国产乱子伦无码精品小说| 国产精品一区二区久久精品无码| 456亚洲人成高清在线| 宅男噜噜噜66国产在线观看| 熟女日韩精品2区| yjizz视频最新网站在线| 视频国产精品丝袜第一页| 一级不卡毛片| 国产一级二级在线观看| 国产激爽爽爽大片在线观看| 狠狠五月天中文字幕| 国产欧美在线视频免费| 91精品啪在线观看国产60岁| 国产日产欧美精品| 欧美日本一区二区三区免费| 国产va在线观看| 99精品热视频这里只有精品7| 亚洲经典在线中文字幕| 色精品视频| 2021国产精品自拍| 亚洲天堂日韩在线| 国产精品夜夜嗨视频免费视频 | 欧美一级特黄aaaaaa在线看片| 欧美日韩资源| 精品视频免费在线| 国产成人综合亚洲欧洲色就色| 在线免费观看a视频| 亚洲一区免费看| 99久久这里只精品麻豆 | 国产精品人莉莉成在线播放| 丁香综合在线| 国产在线日本| 成人免费网站久久久| 欧美一级高清免费a| 国产精品熟女亚洲AV麻豆| 中文字幕无码中文字幕有码在线 | 999福利激情视频 | 久久精品中文无码资源站|