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

動力總成液壓懸置系統固有屬性的拉普拉斯域分析方法*

2016-04-17 06:11:54陳瀟凱
汽車工程 2016年6期
關鍵詞:模態模型系統

宋 康, 陳瀟凱, 林 逸

(1.北京理工大學機械與車輛學院,北京 100081; 2.北京汽車股份有限公司汽車工程研究院,北京 101300)

2016116

動力總成液壓懸置系統固有屬性的拉普拉斯域分析方法*

宋 康1, 陳瀟凱1, 林 逸2

(1.北京理工大學機械與車輛學院,北京 100081; 2.北京汽車股份有限公司汽車工程研究院,北京 101300)

針對液壓懸置隨激振頻率呈非線性變化的動剛度特性,采用了一種拉普拉斯域分析方法計算動力總成液壓懸置系統的固有屬性。首先,建立液壓懸置的集中參數線性流體模型,分析懸置在低頻段和高頻段的穩態動剛度特性,而得到懸置的寬頻動剛度表達式。接著建立懸置系統在曲軸坐標系下的振動模型,在拉普拉斯域內直接應用液壓懸置的動剛度表達式推導系統固有屬性的分析方法。最后采用這一新方法計算了某液壓懸置系統的固有屬性。結果顯示,液壓懸置系統具有復特征值和復模態向量。復特征值同時包含系統的固有頻率和阻尼比信息;根據阻尼比、復模態向量的幅值等信息可對模態類型作出判斷;結合復模態向量的相位信息可確定模態的具體形狀。

動力總成懸置系統;液壓懸置;動剛度;固有屬性;拉普拉斯域

前言

在動力總成懸置系統的設計中,理想的懸置元件應該具有低頻大剛度、大阻尼,高頻小剛度、小阻尼的特性。與彈性體懸置相比,液壓懸置所具有的動剛度特性更加符合該設計要求,因此,液壓懸置在懸置系統的設計中得到了越來越廣泛的應用。

鑒于液壓懸置在特性上的優勢,國內外對液壓懸置開展了廣泛研究,包括對懸置建模、性能分析和對采用液壓懸置的機械系統性能的研究等。文獻[1]中提出了線性時不變集中參數流體-機械模型,該模型在低頻段與試驗數據擬合良好,前人提出的多數模型均可看作該模型在特定條件下的簡化版本。文獻[2]中提出了基于有限測量數據對液壓懸置的性能參數進行預測的方法,分別給出了懸置在低頻段和高頻段的等效機械模型,并在此基礎上提出了液壓懸置的寬頻模型。文獻[3]中建立了慣性通道-解耦盤式液壓懸置的力學模型和數學模型,分析了模型各參數對懸置動特性的影響。文獻[4]和文獻[5]中建立了液壓懸置的線性和非線性集中參數模型并提出了模型參數的識別方法。雖然集中參數模型已能對液壓懸置的頻變剛度特性進行較為精確的預測,但由于該模型具有不同的驅動點剛度和被動點剛度,因此在機械系統的建模中須要對二者加以區別[2,6-7]。文獻[8]中將液壓懸置動剛度的測量值引入車輛或通用機械系統的建模,并根據系統阻尼的不同處理方式總結了特征值問題的多種計算方法。文獻[9]和文獻[10]中在動力總成懸置系統的設計中考慮了液壓懸置內部流體的慣性,增加了振動系統的自由度,建立了擴展系統的時域運動方程,擴展后的模型能同時對動力總成和液壓懸置的固有屬性進行預測。對于動力總成液壓懸置系統,時域方法雖然能夠對系統屬性進行全面預測,但運動方程的建立過程較為繁瑣,并且擴展模型的物理意義不明確;尤其是當擴展系統的維數較大時,計算效率和精度均會下降。

本文中首先建立了液壓懸置的集中參數流體模型,得到懸置動剛度的解析表達式;然后建立了動力總成懸置系統的振動模型,推導得出了系統在拉普拉斯域內的固有屬性分析方法;最后以某液壓懸置系統為例,采用新方法計算得到系統的固有屬性。

1 液壓懸置的集中參數線性模型

以慣性通道-解耦盤式液壓懸置為例,該類型懸置可以簡化為如圖1所示的集中參數流體模型。液壓懸置內部分為2個液室:上液室和下液室。上液室的壓強為p1,體積柔度為C1;下液室的壓強為p2,體積柔度為C2。上、下液室之間通過慣性通道和解耦盤通道相連通。慣性通道的液體慣量為Ii,阻力為Ri,流量為qi;解耦盤通道的液體慣量為Id,阻力為Rd,流量為qd。橡膠主簧的剛度為kr,阻尼為cr,液室的有效作用面積為A。

圖1 液壓懸置集中參數流體模型

液壓懸置內部流體的動力學方程為

(1)

式中:前2個方程為流體的動量方程;后2個方程為流體的連續性方程。

已知液壓懸置的上端與動力總成相連,下端與車身相連。假設懸置上端的位移為x,懸置受到的作用力為F,則F為

(2)

現定義如下:mi=A2Ii,ci=A2Ri,md=A2Id,cd=A2Rd,dxi/dt=qi/A,dxd/dt=qd/A,k1=A2/C1,k2=A2/C2。根據各變量的定義,式(1)與式(2)可改寫為

(3)

方程的改寫表示,由式(1)和式(2)表征的液壓懸置流體模型被轉化為由式(3)表征的機械模型。該等效機械模型在作用效果上與原流體模型相同。

當激勵頻率較低而幅值較大時,解耦盤基本保持關閉狀態,液體主要通過慣性通道進行流動,此時懸置所表現出的性能與慣性通道式液壓懸置相似,懸置動剛度與流體通過慣性通道的能量損失密切相關。當激勵頻率接近慣性通道內液體的固有頻率時,懸置動剛度將達到最大。懸置在低頻段(qd= 0)的動力學方程為

(4)

當激勵頻率較高而幅值較小時,慣性通道內的液體幾乎不再流動,此時液體主要通過解耦盤進行流動。懸置在高頻段(qi= 0)的動力學方程為

(5)

假定系統具有零初始條件(如不作特殊說明,下文在作拉普拉斯變換時均假定初始條件為0)并考慮到k1遠大于k2,對式(4)和式(5)進行拉普拉斯變換得到懸置的動剛度表達式為

(6)

對于式(6),在低頻段內,m0=mi,c0=ci;在高頻段內,m0=md,c0=cd。式(6)對應液壓懸置的寬頻模型,研究結果表明,該模型的有效頻率范圍可達300Hz[2]。

需要說明,上述給出的液壓懸置集中參數流體模型或等效機械模型屬于穩態模型,該模型只能針對低頻激勵或高頻激勵作穩態分析,但無法對頻率的過渡段給出可信的結果,這是由于該頻段涉及解耦盤工作狀態的選擇機制。另外,建模過程中還忽略了諸如真空效應、湍流等典型的非線性效應[11]。

2 動力總成懸置系統的固有屬性分析方法

對于采用彈性體懸置的動力總成懸置系統,在小阻尼條件下,系統的固有屬性可根據質量(慣量)特性和剛度特性進行確定,計算得到的特征值為實數,特征向量為實向量。如果考慮阻尼對系統振動特性的影響,則確定系統的固有屬性時還要用到阻尼特性,此時計算得到的特征值一般為復數,特征向量為復向量。

對于彈性體懸置,動剛度和滯后角基本隨頻率呈線性變化,這就說明懸置的剛度kr和阻尼cr基本可以是頻率恒定的。彈性體懸置動剛度的表達式為

(7)

對于液壓懸置,隨著激勵頻率的增加,懸置的動剛度和滯后角均呈現出典型的非線性變化關系,因此,對于采用液壓懸置的動力總成懸置系統,其固有屬性的分析將更加復雜。針對該問題,本節提出了懸置系統在拉普拉斯域內的固有屬性分析方法。

2.1 系統建模

在對動力總成懸置系統建模時,將動力總成簡化為剛體;彈性體懸置與液壓懸置均簡化為3個相互正交的彈簧,各彈簧剛度以動剛度的形式給出。動力總成懸置系統的振動模型如圖2所示。

圖2 動力總成懸置系統的振動模型

以曲軸坐標系{XC,YC,ZC}作為建立系統運動方程和計算固有屬性的參考坐標系。曲軸坐標系的原點位于動力總成質心處,XC軸沿曲軸軸線指向汽車前方,ZC軸沿直列氣缸軸線或V型氣缸角平分線指向上方,YC軸根據右手定則確定。同時,在每個懸置處建立局部坐標系,第i個懸置處的局部坐標系定義為{XLi,YLi,ZLi},懸置的動剛度在局部坐標系中表示。

2.2 拉普拉斯域內的固有屬性分析方法

在拉普拉斯域內,液壓懸置在流體起作用方向(通常定義為垂向)的動剛度與頻率呈典型的非線性關系,而在另外兩個方向的動剛度則近似認為與彈性體懸置相同。在局部坐標系內,懸置的動剛度可表示為

(8)

5)關聯系數:式中:εi(k)為X0與Xi在第k個指標的關聯系數;|X0(k)-Xi(k)|為與在第k個指標的絕對值;min min|X0(k)-Xi(k)|為二級最小絕對差;max max|X0(k)-Xi(k)|為二級最大絕對差;ρ為分辨系數,一般取值0.5。

重新觀察式(8),對于液壓懸置的垂向動剛度表達式,如果將公式中的參數m0與c0設置為0,則該式在形式上就與彈性體懸置的動剛度表達式相同。因此,通過對液壓懸置剛度矩陣中的參數作不同設置,可以得到彈性體懸置與液壓懸置剛度矩陣的通用表達式。根據液壓懸置的物理結構,將參數m0與c0設置為0是忽略懸置內部流體的作用,這種情況相當于去除液壓懸置的液壓油,只有橡膠主簧發揮作用。

根據旋量理論,第i個懸置局部坐標系原點的位移dLi與動力總成質心的位移d滿足如下關系式:

(9)

式中:RC,Li為由第i個局部坐標系至曲軸坐標系的旋轉矩陣;ri為第i個局部坐標系在曲軸坐標系中的位置向量;ad(ri)為向量ri的伴隨變換形式。

對式(9)等號左右兩側同時進行拉普拉斯變換L可得

(10)

根據第i個懸置在局部坐標系中的位移dLi及由此產生的反作用力FLi計算得到懸置動剛度KLi(s)的表達式為

(11)

其中動剛度KLi(s)表示為

根據旋量理論,作用在動力總成質心且與懸置作用力FLi等效的作用力Fi表示為

Fi=AdgC,LiFLi

(12)

對所有力旋量Fi進行求和計算,所得結果即為動力總成質心的位移向量d引起的作用力旋量F:

(13)

對式(13)等號左右兩側同時進行拉普拉斯變換可得

(14)

將式(10)~式(12)代入式(14)可得

(15)

根據式(15)可得系統動剛度矩陣K(s)的表達式為

(16)

動力總成懸置系統的運動方程可寫為

(17)

式中:M為系統的廣義質量矩陣;FE為系統受到的外部激勵。

對式(17)等號左右兩側同時進行拉普拉斯變換并代入式(15)和式(16)可得

s2ML[d]+K(s)L[d]=L[FE]

(18)

Z(s)=s2M+K(s)

系統傳遞函數T(s)的表達式為

(19)

式中adj[Z(s)]表示Z(s)的伴隨矩陣。

根據振動模態理論,式(19)的分母,即Z(s)的行列式,稱為系統的特征方程。特征方程的根,即為系統的特征值,每個特征值都同時包含系統的固有頻率與模態阻尼比信息[12]。

對于系統的第r個特征值λr,記與之相應的特征向量為ur。特征值λr與特征向量ur滿足如下方程:

(20)

根據式(20)計算得到的特征向量只含有6個元素,對應動力總成的6個自由度。對于液壓懸置的等效機械模型,當系統發生固有振動時,流體集中質量元素也具有相應的位移xM,i(xi或xd)。對于系統的第r階固有振動,根據式(6)計算得到位移xM,i:

(21)

式中:uLi,r為第i個懸置在系統發生第r階固有振動時的位移向量,該向量在局部坐標系中表示,可參照式(10)坐標系轉換方法進行計算;uLi,r(3)為向量的第3個元素。

對于采用了N個液壓懸置的動力總成懸置系統,完整的系統特征向量Ur可寫為

(22)

至此,根據式(19)和式(22)可分別計算得出系統的全部特征值和完整的特征向量。由于系統的傳遞函數在拉普拉斯域中表示,因此推導過程中可直接采用式(6)給出的液壓懸置動剛度表達式。與時域方法相比,傳遞函數矩陣仍然為6維矩陣,但是卻可以一次性計算得到6+N個特征值和相應的特征向量,無須對振動系統的自由度進行擴展。而且,式(21)具有明確的物理意義,反映了當系統發生固有振動時,流體集中質量元素的位移與動力總成位移之間的傳遞關系。

3 應用算例

以某2.3L直列四缸四沖程動力總成的懸置系統為例,說明拉普拉斯域分析方法在系統固有屬性計算中的應用。已知動力總成的質量為216kg,慣量參數在曲軸坐標系中測量,結果如表1所示。

表1 動力總成慣量參數 kg·m2

系統共采用3處懸置,分別為安裝于變速器一側的左懸置、安裝于發動機一側的右懸置和安裝于發動機后方的后懸置。其中,左懸置與后懸置為彈性體懸置,右懸置為慣性通道式液壓懸置。系統中各懸置的安裝位置、安裝方向和彈性體懸置的剛度值分別如表2~表4所示。懸置的安裝位置在曲軸坐標系中描述;安裝方向以局部坐標系與曲軸坐標系各坐標軸之間的夾角表示;懸置的剛度和阻尼在局部坐標系中描述。

表2 懸置安裝位置 mm

表3 懸置安裝方向 (°)

表4 懸置剛度 N·mm-1

液壓懸置的結構參數為:橡膠主簧剛度kr= 168N/mm,主簧阻尼cr= 100N·s/m,流體等效集中質量mi= 50.3kg,等效剛度ki= 250N/mm,等效阻尼ci= 2582N·s/m。另外,所有彈性體懸置的阻尼均設置為100N·s/m。

根據第2節給出的動力總成液壓懸置系統在拉普拉斯域內的固有屬性分析方法,計算系統的特征值與特征向量(模態向量)。對計算得到的特征值作進一步處理而得到系統的固有頻率和模態阻尼比,結果如表5所示;由于模態向量以各元素(復數)比值的形式給出,因此設定向量中第1個元素(對應動力總成的側向平動自由度)的幅值為1,相位角為0,參考該設定計算得出向量中的其他各元素。系統復模態向量的幅值和相位信息分別如表6和表7所示。在表6中,Tx,Ty和Tz分別表示模態向量幅值在XC軸,YC軸和ZC軸平動自由度上的分量;Rx,Ry和Rz分別表示向量幅值在相應各軸轉動自由度上的分量;Tm表示向量幅值在流體等效集中質量平動自由度上的分量。

表5 系統各階固有頻率與阻尼比

表6 系統復模態向量的幅值

表7 系統復模態向量的相位 (°)

該動力總成懸置系統具有7階固有頻率。其中,6階固有頻率對應動力總成剛體的共振頻率,1階固有頻率對應液壓懸置慣性通道內流體的共振頻率。

根據表5中的結果可知,第1階模態的阻尼比為40.53%,遠大于其他各階模態的阻尼比;同時,根據表6中的結果可知,在第1階模態向量中,流體等效集中質量的相對幅值Tm為6.61,是向量中的最大分量。因此,根據上述分析可以判定,系統的第1階固有頻率對應液壓懸置慣性通道內流體的共振頻率;由于懸置內部的等效阻尼非常大,因此系統的阻尼比處于相對較高的水平。

根據表6中的幅值信息可知,系統的第4階、第6階和第7階模態向量中均存在幅值占比較大的元素,因此判定這3階模態分別對應系統的俯仰模態、橫擺模態和側傾模態。其中,系統的第4階模態尤為重要,這是因為該階模態與系統在轉矩作用下的隔振性能設計密切相關。

對于系統的第2階、第3階和第5階模態,無法根據幅值信息對模態類型作出直觀判斷。例如,對于第2階模態向量,雖然最大元素對應Rx方向,但與第4階模態相比,Rx元素在模態向量中所占的比重顯然要小一些;對于第3階和第5階模態,模態向量中存在多個幅值相差不大的元素,因此不易作出直觀判斷。

為了對模態的類型進行判定,現提出2條參考標準:第1條標準是模態能量百分比。由于動力總成的質量與慣量在大小上基本相差一個數量級,因此,模態向量中的平動分量對模態能量百分比的影響更大。參考這一標準,判定第2階模態對應系統的縱向平動模態,第3階模態對應系統的側向平動模態,第5階模態對應系統的垂向跳動模態。

第2條標準是液壓懸置集中流體質量的運動與動力總成各向運動之間的耦合程度。由于液壓懸置安裝于動力總成右側,因此,集中流體質量的運動將與動力總成的垂向跳動和側傾轉動產生較強的耦合,這一點可以根據第1階模態向量和第7階模態向量進行驗證。前述分析已經判定第1階頻率對應液壓懸置慣性通道內流體的共振頻率,根據表6中的信息可知,第1階模態向量中除了最大的Tm分量外,Tz分量和Ry分量的幅值也較大;對于第7階模態,已經判定該階模態對應系統的側傾模態,模態向量中除了最大的Ry分量,Tm分量的幅值也較大。第1階模態和第7階模態的分析結果符合第2條標準。因此,對于系統的垂向跳動模態,模態向量中的Tz分量和Tm分量一定同時具有較大的幅值,據此判定第5階模態對應系統的垂向跳動模態。

根據表6中的幅值信息雖然可以對模態的類型作出判斷,但是還需要根據表7中的相位信息確定模態的具體形狀。相位信息對系統第1階、第5階和第7階模態形狀的判斷尤為重要。對于集中流體質量的運動和動力總成的垂向跳動,雖然由第1階和第5階模態向量中的幅值信息可以看出二者的耦合效應非常明顯,但根據表7中的相位信息可知,二者在第1階模態中的相位差為15.54°,而在第5階模態中的相位差為133.33°,因此,二者在系統發生第1階固有振動時傾向于同相運動,而在第5階固有振動時傾向于反相運動。同樣,對于集中流體質量的運動和動力總成的側傾轉動,二者在第1階固有振動時傾向于反相運動,而在第7階固有振動時傾向于同相運動。

4 結論

(1)液壓懸置的集中參數流體模型或等效機械模型可較為準確地反映懸置的頻變動剛度特性,能夠同時用于懸置在低頻段與高頻段的穩態性能分析;但由于該模型為穩態模型,因此無法模擬懸置的解耦盤選擇機制等非線性行為。

(2)拉普拉斯域分析方法可直接將液壓懸置的動剛度表達式用于系統傳遞函數的推導;在不增加系統自由度的前提下,該方法可以一次性計算出系統的全部特征值;推導過程同時涉及系統的總體剛度特性和液壓懸置的結構特性,因而物理意義明確。

(3)動力總成液壓懸置系統具有復特征值及復模態向量。復特征值同時包含系統的固有頻率及阻尼比信息;根據阻尼比、復模態向量的幅值等信息可以對模態類型作出判斷;結合復模態向量的相位信息可以確定模態的具體形狀。

[1]SINGHR,KIMG,RAVINDRAPV.Linearanalysisofautomotivehydro-mechanicalmountwithemphasisondecouplercharacteristics[J].JournalofSoundandVibration, 1992, 158(2): 219-243.

[2]HES,SINGHR.Estimationofamplitudeandfrequencydependentparametersofhydraulicenginemountgivenlimiteddynamicstiffnessmeasurements[J].NoiseControlEngineeringJournal, 2005, 53(6): 271-285.

[3] 馮振東,王立公,宋傳學.動力總成液力懸置的動特性仿真[J].汽車工程,1994,16(4):193-198.

[4] 呂振華,梁偉,上官文斌.汽車發動機液阻懸置動特性仿真與實驗分析[J].汽車工程,2002,24(2):105-111.

[5] 上官文斌,呂振華.液阻型橡膠隔振器非線性特性仿真分析[J].振動工程學報,2003,16(4):393-398.

[6]LEEJH,BAEMS,KIMKJ.Limitationsofmechanicalmodelwithlumpedmassinrepresentingdynamiccharacteristicsofhydraulicmount[C].SAEPaper2003-01-1466.

[7]LEEJH,SINGHR.Criticalanalysisofanalogousmechanicalmodelsusedtodescribehydraulicenginemounts[J].JournalofSoundandVibration, 2008, 311(3-5): 1457-1464.

[8]JEONGT,SINGHR.Inclusionofmeasuredfrequencyandamplitude-dependentmountpropertiesinvehicleormachinerymodels[J].JournalofSoundandVibration, 2001, 245(3): 385-415.

[9]PARKJY,SINGHR.Roleofspectrallyvaryingmountpropertiesininfluencingcouplingbetweenpowertrainmotionsundertorqueexcitation[J].JournalofSoundandVibration, 2010, 329(14): 2895-2914.

[10]HUJF,CHENWW,HUANGH.Decouplinganalysisforapowertrainmountingsystemwithacombinationofhydraulicmounts[J].ChineseJournalofMechanicalEngineering, 2013, 26(4): 737-745.

[11]ADIGUNAH,TIWARIM,SINGHR,etal.Transientresponseofahydraulicenginemount[J].NoiseControlEngineeringJournal, 2003, 268(2): 217-248.

[12] 沃德·海倫,斯蒂芬·拉門茲,波爾·薩斯.模態分析理論與試驗[M].北京:北京理工大學出版社,2001:3-26.

Analysis Method of Natural Properties of Powertrain Hydraulic Mounting System in Laplace Domain

Song Kang1, Chen Xiaokai1& Lin Yi2

1.SchoolofMechanicalEngineering,BeijingInstituteofTechnology,Beijing100081; 2.BeijingAutomotiveTechnologyCenter,BAICMotorCorporationLtd.,Beijing101300

In view of that the dynamic stiffness of hydraulic mount nonlinearly changes with excitation frequency, a Laplace domain analysis method is adopted to calculate the natural properties of powertrain hydraulic mounting system. Firstly a linear fluid model with lumped parameters for hydraulic mount is built, the steady state dynamic stiffness characteristics of powertrain mount in both low and high frequency bands are analyzed, and the dynamic stiffness expression for powertrain mount in a wide frequency range is obtained. Then a vibration model for powertrain mount is created in crankshaft coordinates and an analysis method of the natural properties of system is derived by directly applying dynamic stiffness expression of hydraulic mount in Laplace domain. Finally this new method is used to calculate the natural properties of a hydraulic mounting system. The results indicate that hydraulic mounting system has complex eigenvalues and complex modal vectors, in which eigenvalues contain information of the natural frequency and damping ratio of system. The type of mode can be judged according to damping ratio and the magnitude of complex modal vector, while the specific mode shapes can be determined based on phase information of complex modal vectors.

powertrain mounting system; hydraulic mount; dynamic stiffness; natural properties; Laplace domain

*國家自然科學基金(51275040)資助。

原稿收到日期為2014年12月25日。

猜你喜歡
模態模型系統
一半模型
Smartflower POP 一體式光伏系統
工業設計(2022年8期)2022-09-09 07:43:20
WJ-700無人機系統
ZC系列無人機遙感系統
北京測繪(2020年12期)2020-12-29 01:33:58
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
連通與提升系統的最后一塊拼圖 Audiolab 傲立 M-DAC mini
3D打印中的模型分割與打包
國內多模態教學研究回顧與展望
基于HHT和Prony算法的電力系統低頻振蕩模態識別
主站蜘蛛池模板: 日韩在线影院| 国产91小视频在线观看| av一区二区三区高清久久| 91青青草视频| 久久久久青草线综合超碰| 国产成人精品第一区二区| 婷五月综合| 四虎永久在线视频| 国产成人一级| 日韩123欧美字幕| 免费99精品国产自在现线| 欧美日韩中文国产| 日韩二区三区无| 亚洲精品少妇熟女| 无码国产伊人| 青青青国产视频| 欧美成人日韩| 亚洲成aⅴ人在线观看| 国产精品福利一区二区久久| 午夜国产精品视频黄| 久久一级电影| 91精品啪在线观看国产91九色| 很黄的网站在线观看| 亚洲高清无在码在线无弹窗| 国产精品欧美亚洲韩国日本不卡| 日本成人一区| 精品一区二区三区水蜜桃| 亚洲VA中文字幕| 不卡无码网| 日本人妻一区二区三区不卡影院 | 亚洲av无码专区久久蜜芽| 91在线一9|永久视频在线| 国产精品太粉嫩高中在线观看 | 欧美激情视频二区三区| 国产国产人成免费视频77777 | 日韩国产高清无码| 亚洲高清在线天堂精品| 日韩欧美国产精品| www成人国产在线观看网站| 亚洲精品成人福利在线电影| 超薄丝袜足j国产在线视频| 中文无码精品a∨在线观看| 69av免费视频| 国产在线精品99一区不卡| 国产成人综合亚洲欧美在| 日本在线免费网站| 成人a免费α片在线视频网站| 久无码久无码av无码| 欧美精品综合视频一区二区| 亚洲 欧美 偷自乱 图片 | 91免费在线看| 久久精品国产电影| 免费无码AV片在线观看中文| 亚洲无码91视频| 2021国产乱人伦在线播放| 五月激情婷婷综合| 午夜免费视频网站| 精品少妇人妻av无码久久| 亚洲愉拍一区二区精品| 国产在线小视频| 午夜小视频在线| 国产精品第页| 22sihu国产精品视频影视资讯| www.99精品视频在线播放| 国产精品永久久久久| 91精品伊人久久大香线蕉| 无遮挡国产高潮视频免费观看| 国产精品99r8在线观看| 日韩中文无码av超清| 无码精品一区二区久久久| 一本无码在线观看| 日韩在线2020专区| 亚洲国产系列| 九九九精品视频| 亚洲国产天堂久久综合| 亚洲AV无码久久天堂| 在线欧美一区| 精品国产Av电影无码久久久| 日韩中文欧美| 在线观看av永久| 久久77777| 欧美精品三级在线|