亚洲免费av电影一区二区三区,日韩爱爱视频,51精品视频一区二区三区,91视频爱爱,日韩欧美在线播放视频,中文字幕少妇AV,亚洲电影中文字幕,久久久久亚洲av成人网址,久久综合视频网站,国产在线不卡免费播放

        ?

        波浪地震耦合作用下懸浮隧道動力響應(yīng)分析

        2021-03-02 09:59:04張玉龍潘少康
        工程力學(xué) 2021年2期
        關(guān)鍵詞:管體波浪錨索

        羅 剛,張玉龍,潘少康,任 毅

        (1.長安大學(xué)公路學(xué)院,西安710064;2.廣西新發(fā)展交通集團(tuán)有限公司,南寧530029;3.紹興文理學(xué)院土木工程學(xué)院,紹興312000)

        懸浮隧道又稱阿基米德橋是挪威學(xué)者提出的一種跨越長大水域的新型結(jié)構(gòu)形式[1],相比于其它結(jié)構(gòu)由于其懸浮在水中一定深度使其具有受外界環(huán)境干擾小、造價低、綠色環(huán)保等優(yōu)點[2]。目前懸浮隧道尚處于研究階段,包括挪威、意大利、中國、韓國、日本等國學(xué)者都對其開展了大量的研究[3-8]。

        波浪荷載是懸浮隧道服役期間受到的長期荷載。在波浪荷載的計算中,Morison 方程被認(rèn)為是一種簡便的計算方法,項貽強(qiáng)等[9]通過Stokes波浪理論和Morison 方程研究了線性波浪理論計算懸浮隧道波浪力的適用范圍。葛斐等[10-11]通過Airy 線性波浪理論和Morison 方程對懸浮隧道錨索的在波浪荷載下的動力響應(yīng)和海流作用下的渦激振動進(jìn)行了研究。麥繼婷等[12]通過Morison 方程分析了懸浮隧道放置深度、海流速度、波浪周期等因素對水平波浪力的影響。對于更為復(fù)雜的非線性波浪力問題則可以采用勢流理論進(jìn)行分析,Paik 等[13]采用邊界元法將懸浮隧道簡化為三維梁單元,基于繞射理論求解無粘不可壓縮的Navier-Stokes方程對懸浮隧道受到的波浪力進(jìn)行了研究。Wu 等[14]基于勢流理論針對懸浮隧道受到的非線性波浪力進(jìn)行了研究。除了荷載的計算外,波浪荷載下懸浮隧道的響應(yīng)問題也是研究的重點,Kunisu[15]采用邊界元法對懸浮隧道在波浪荷載下的動力響應(yīng)進(jìn)行了計算,分析了隧道尺寸和形狀對其響應(yīng)的影響。Remseth 等[16]提出了一種計算波浪荷載下懸浮隧道隨機(jī)動力響應(yīng)的方法,分析了波向、阻尼系數(shù)、剛度和自振頻率的影響。Faggiano等[17]對波浪力作用下懸浮隧道結(jié)構(gòu)的響應(yīng)進(jìn)行詳細(xì)分析,討論的懸浮隧道不同方向的受力情況以及隨深度的變化規(guī)律。Dai等[18]以帶錨索的浮筒式懸浮隧道為研究對象,分析了長期波浪荷載作用下懸浮隧道的響應(yīng)問題。Jin等[19]采用室內(nèi)試驗和數(shù)值模擬相結(jié)合的方式分析了波浪作用下黏滯阻尼系數(shù)、浮重比、錨索剛度對懸浮隧道動力響應(yīng)的影響。

        懸浮隧道的修建場地往往在有地震帶的海峽上,這給懸浮隧道帶來極大的威脅。在此方面Morita 等[20]以格林函數(shù)為基礎(chǔ)將懸浮隧道錨索簡化為支撐彈簧,進(jìn)行了豎向地震作用下懸浮隧道的響應(yīng)分析。Lee等[21]在考慮海水的可壓縮性、彈性海床對壓力波吸收的條件下,研究了海水對SFT系統(tǒng)地震響應(yīng)的主要影響。數(shù)值模擬是研究地震的重要手段,Leira[22]采用大質(zhì)量法對多點激勵下SFT 的動力響應(yīng)進(jìn)行了數(shù)值模擬,研究了地震波通道效應(yīng)對懸浮隧道的動力特性的影響。Fogazzi等[23]采用改進(jìn)的桿單元來模擬錨索,梁單元模擬管體建立懸浮隧道的有限元模型,進(jìn)行了地震激勵下的動力響應(yīng)分析。上述研究多是針對單一荷載進(jìn)行討論,然而工程實際中荷載往往都不是單一出現(xiàn)的,因此有一些學(xué)者開展了多荷載耦合下懸浮隧道的響應(yīng)研究。Jin 和Kim[24]在波浪荷載的基礎(chǔ)上考慮了移動荷載,研究了兩者耦合作用下懸浮隧道的特性。DiPilato等[25]采用數(shù)值模擬的方式對波浪、地震耦合作用下懸浮隧道的動力特性進(jìn)行了研究。Wu 等[26]基于三角級數(shù)法得到了隨機(jī)地震荷載,建立了地震作用下的水動力荷載計算公式,對懸浮隧道錨索在地震和流體共同作用下的動力特性進(jìn)行了研究。此外,一些學(xué)者也對其他海洋結(jié)構(gòu)進(jìn)行了多荷載耦合作用的響應(yīng)分析[27-28]。

        對懸浮隧道進(jìn)行合理的簡化,建立合適的理論模型也是其動力分析的重要部分。孫勝男等[29-30]將懸浮隧道錨索簡化為簡支梁,分析了錨索在軸向激勵和隨機(jī)激勵下的動力響應(yīng)。惠磊等[31]等建立了沖擊作用下懸浮隧道的圓柱殼模型,采用等效質(zhì)量法進(jìn)行求解,分析了碰撞中的能量轉(zhuǎn)化和沖擊點的最大位移。陳健云等[32-33]通過將管體簡化為質(zhì)量塊,將錨索簡化為歐拉梁建立了懸浮隧道耦合振動模型,考慮了水流作用下懸浮隧道的耦合振動響應(yīng)。Lin 等[34]通過將管體和錨索簡化為耦合連接的歐拉梁系統(tǒng),建立了懸浮隧道管體-錨索耦合振動模型,研究了移動荷載下懸浮隧道的動力響應(yīng)。總體而言,目前關(guān)于懸浮隧道研究主要集中在單一荷載、特定構(gòu)件的研究,對于荷載耦合和結(jié)構(gòu)耦合的研究較少,不能較為真實反映懸浮隧道的實際工況。

        為了綜合考慮波浪地震荷載以及懸浮隧道結(jié)構(gòu)的耦合效應(yīng),本文在前人研究的基礎(chǔ)上,通過Stokes波浪方程和三角級數(shù)法計算波浪荷載和地震荷載,基于D’Alembert 原理建立波浪地震作用下管體-錨索耦合振動模型,通過四階Runge-Kutta 法數(shù)值積分求解微分方程組,對地震、波浪的荷載參數(shù)和系統(tǒng)響應(yīng)進(jìn)行分析和討論。

        1 荷載模型

        1.1 波浪荷載

        采用微幅波理論的一階近似不能滿足波浪荷載計算精度要求,應(yīng)當(dāng)考慮非線自由表面對波浪荷載的影響。Stokes波浪理論是計算非線性波浪荷載的常用方法,在保證計算精度的情況下,本文選擇Stokes三階波浪方程來進(jìn)行計算,波面方程、速度勢函數(shù)和彌散關(guān)系為:

        在已知波高H、周期T以及水低深度d的情況下,可以先求波長L和a,進(jìn)一步則可解出Stokes三階波浪方程的全部方程。將波浪引起的水質(zhì)點速度式(2)代入Morison 方程即可計算出作用在懸浮隧道上的波浪荷載,在后文中將進(jìn)行這部分計算。

        1.2 地震荷載

        為了對懸浮隧道進(jìn)行地震響應(yīng)時域分析,必須建立合適的地震加速度時程曲線。一種可行的方法是通過三角級數(shù)法計算出平穩(wěn)過程條件下地震加速度時程曲線[35]:

        式中:S(ω)為設(shè)計地震功率譜;Δω為圓頻率增量,Δω=2π/Td,Td是地震的持續(xù)時間;ωk=kΔω;εk為隨機(jī)相位改變量,其變化范圍為0~2π。

        考慮到隨機(jī)地震激勵的非平穩(wěn)特性,根據(jù)功率譜包絡(luò)理論中的非均勻調(diào)制函數(shù)對式(3)進(jìn)行修正[36]:

        式中:t0和tn分別為控制平穩(wěn)段的首末時刻;c為衰減常數(shù)。

        對于地震譜的選擇具有較大任意性,結(jié)合實際的場地條件和結(jié)構(gòu)特性采用設(shè)計地震加速度響應(yīng)譜,能夠使懸浮隧道的抗震分析更加的精細(xì)化和合理。目前中國罕有海洋結(jié)構(gòu)相關(guān)的抗震設(shè)計標(biāo)準(zhǔn)。因此,本文參考國際標(biāo)準(zhǔn)化組織(ISO)制定的海洋結(jié)構(gòu)抗震設(shè)計標(biāo)準(zhǔn)[37],地震設(shè)計反應(yīng)譜需要根據(jù)場地分類、峰值加速度和阻尼比三個參數(shù)進(jìn)行選擇,具體如圖1所示。

        圖1 地震設(shè)計加速度響應(yīng)譜Fig.1 Seismic design acceleration response spectrum

        通過設(shè)計地震響應(yīng)譜則地震的功率譜可以近似通過下式得到[38]:

        式中:ζ 為結(jié)構(gòu)的阻尼比,此處取0.05;Sa(ω)為與Sa(T)對應(yīng)的圓頻率響應(yīng)譜,可由對后者進(jìn)行T=2π/ω變量替換得到;P為超越概率常數(shù),在此取0.9。

        由于式(5)的轉(zhuǎn)換存在誤差,需要對人工地震功率譜進(jìn)行修正。對于每一個頻率點處的人工地震功率譜,以計算反應(yīng)譜和目標(biāo)反應(yīng)譜的比值作為驗證精度值,當(dāng)精度值不滿足擬合精度值時,按照比值進(jìn)行功率譜值的調(diào)整,調(diào)整后再次生成人工地震動,并用新人工地震動計算相應(yīng)的反應(yīng)譜,再與目標(biāo)反應(yīng)譜比較,進(jìn)行精度值檢驗和調(diào)整。如此反復(fù),直至達(dá)到精度要求。

        通過上述方法,可以建立與目標(biāo)場地和結(jié)構(gòu)相匹配的加速度時程曲線,則作用在結(jié)構(gòu)上的地震荷載可用下式來表示:

        式中:Fg為作用在結(jié)構(gòu)上的地震荷載;m為結(jié)構(gòu)的質(zhì)量。

        1.3 流體作用

        如果流體在某一點處的速度和加速度已知,則可以通過Morison 方程即式(7)計算出流體的作用力,方程的第一項是流體作用引起的阻尼力,第二項是流體作用引起的慣性力。

        式中:ρw海水的密度;D為懸浮隧道的外徑;v為水質(zhì)點速度;w為管體振動的位移函數(shù);CD和CM分別為流體作用的阻尼力系數(shù)和附加質(zhì)量系數(shù),可通過流體試驗測定,在缺乏大量試驗的條件下CD可近似取1.0,CM可近似取2.0[9]。

        已有研究表明:在相同荷載作用下,懸浮隧道水平方向的響應(yīng)遠(yuǎn)小于豎直方向[39],因此本文僅考慮懸浮隧道的豎向受力及響應(yīng)。波浪和地震荷載會通過流體間接的作用在懸浮隧道系統(tǒng)上,結(jié)構(gòu)自身的運動也會導(dǎo)致流體作用,此處將同時通過Morison 方程計算三者的流體作用。在計算中錨索受到的流體作用與隧道管體形式相同,但考慮到錨索的直徑較小且具有較大的安裝水深,在此假設(shè)錨索不受到波浪力的作用[11]。

        式中:u為錨索的位移;di為錨索的直徑;xg為地震引起的水質(zhì)點位移。

        2 隧道模型及方程求解

        目前在懸浮隧道理論計算方面,主要有彈性地基梁模型(Beam on elastic foundation)、錨索振動模型(Tether vibration model)和管體-錨索耦合振動模型(Tube-tether coupled model)三種[40-41]。其中,彈性地基梁模型將懸浮隧道錨索簡化為彈性支撐,將管體簡化為彈性支承梁,能夠?qū)荏w的動力響應(yīng)進(jìn)行分析但忽略了錨索的振動。錨索振動模型將懸浮隧道錨索簡化為傾斜的梁,將管體簡化為質(zhì)量塊或作用在錨索上的簡諧激勵,能夠?qū)﹀^索的動力響應(yīng)進(jìn)行分析,但忽略了管體的振動。管體-錨索耦合振動模型同時考慮了錨索和管體的振動以及耦合效應(yīng)能夠同時對管體和錨索進(jìn)行動力響應(yīng)分析,因此本文將選用該模型進(jìn)行計算。

        2.1 錨索振動

        如圖2(a),在懸浮隧道的兩側(cè),假設(shè)等間距h布置著Ns個相同的錨索,管體在浮力作用下使錨索處于繃緊的狀態(tài),此時錨索具有一定的初張力T0。管體發(fā)生w的位移后,在與錨索連接的xi處與其發(fā)生位移協(xié)調(diào),同時錨索自身也會產(chǎn)生的法向振動,因此,除初張力外錨索還承受附加的張力ΔT。

        此外,由于錨索自身的重力和柔度較大,往往表現(xiàn)出一定的垂度,使其受力特性發(fā)生改變。為了考慮錨索的垂度效應(yīng),采用等效剛度法來計算錨索的彈性模量[42]。

        2.2 系統(tǒng)微分方程及求解

        懸浮隧道管體受到波浪力、流體作用力以及錨索的張力,通過錨索的張力項隧道管體與錨索之間相互耦合,在錨索振動方程的基礎(chǔ)上可以寫出懸浮隧道管體-錨索系統(tǒng)的耦合振動方程,如式(11)。

        式中:下標(biāo)i代表第i根錨索,mi為單位長度錨索的質(zhì)量;ui為錨索的位移;ci為錨索的粘滯阻尼系數(shù);Fdi為第i根錨索受到的流體作用力。

        圖2 懸浮隧道耦合振動系統(tǒng)Fig.2 Coupled vibration system of SFT

        式中:mtb為單位長度管體的質(zhì)量;ctb為管體的粘滯阻尼系數(shù);EtbItb為管體的彎曲剛度;δ 為狄拉克函數(shù)。

        上述微分方程組可采用Galerkin 法即式(12)進(jìn)行化簡,將偏微分方程組轉(zhuǎn)換為常微分方程組:由于在錨索的振動中第一階振型占主導(dǎo)地位[43],在此取錨索的第一階模態(tài)進(jìn)行簡化計算,化簡后的結(jié)果為:

        上述的時變耦合微分方程組無法求得解析解,對于懸浮隧道管體微分方程而言,計算n階振型需要聯(lián)立求解n個微分方程。對于懸浮隧道錨索而言,錨索數(shù)量Ns對應(yīng)著需要求解Ns個微分方程。本文計算對耦合微分方程系統(tǒng)進(jìn)行降階處理,然后采用四階龍格庫塔法數(shù)值積求解,具體步驟如圖3所示。

        圖3 系統(tǒng)微分方程系統(tǒng)求解步驟Fig.3 Stepsof solving system differential equation

        參考國內(nèi)外待建懸浮隧道的參數(shù),在下文的計算中所用參數(shù)如表1。

        3 計算與分析

        3.1 不同地震譜

        為了研究不同場地類別和地震作用下懸浮隧道的動力響應(yīng),在此根據(jù)已經(jīng)建立的設(shè)計地震譜選擇了三條與之匹配的地震記錄進(jìn)行對比分析,如表2所示,同時將四種地震的加速的時程曲線如圖4(a)所示。

        表1 懸浮隧道系統(tǒng)基本參數(shù)Table 1 Primary parametersof SFTsystem

        表 2天然地震波參數(shù)Table 2 Parameters of natural earthquakes

        由圖4(a)可以看出,基于三角級數(shù)法生成的人工地震波,在時域上具有較大的分布密度和加速度峰值,主要地震作用時間在前40 s。而基于設(shè)計譜在太平洋地震工程研究中心(Peer)上選擇的地震波具有較大的差異性,主要作用時間及加速度幅值均不相同,但能為懸浮隧道的抗震性能進(jìn)行多方面的檢驗。圖4(b)為本文中錨索位移的計算結(jié)果與Wu 等[26]和Su[41]的結(jié)果進(jìn)行比較,在計算的過程選擇了相同的地震荷載參數(shù)和錨索結(jié)構(gòu)參數(shù)。后者采用的計算模型為錨索振動模型,該模型將懸浮隧道管體的作用力簡化為作用在錨索上的簡諧荷載,主要考慮錨索的振動作用。由于采用了隨機(jī)生成的地震波,計算出來的結(jié)果必然存在一定的差異。盡管本文所用懸浮隧道模型與他們不同,但計算結(jié)果較為相似。值得注意的本文在此對比的是3#錨索的位移,2#錨索因為處于管體的跨中位置,同時受到參數(shù)共振的影響產(chǎn)生了較大位移,這在非耦合模型中是沒有體現(xiàn)的。

        通過多點虛擬激勵法將人工地震和自然地震荷載施加在懸浮隧道系統(tǒng)上,同時考慮波浪荷載計算后得到圖4(c)~圖4(d)動力響應(yīng)結(jié)果。由圖4(c)可知,在人工地震作用下懸浮隧道管體的峰值位移為0.30 m,出現(xiàn)在11.7 s。人工地震波的持續(xù)時間約為40 s,從圖中可以看出前40 s,管體震蕩較為劇烈。另外從圖中也能看出2#錨索位移最大,為2.15 m,出現(xiàn)在9.9 s,其余兩對錨索的最大位移分別為0.57和0.27 m。

        由圖4(d)可以發(fā)現(xiàn),不同地震荷載輸入下結(jié)構(gòu)的響應(yīng)具有較大的差異性,Helena Montana-01地震的特點是峰值加速度大、持續(xù)時間短,因此前20 s內(nèi)管體的振動最為劇烈,但所造成的位移不大。Imperial Valley-02地震與前者相似,但地震持續(xù)時間較長,峰值加速度也更大,其管體最大位移可達(dá)0.41 m。Humbolt Bay 地震的水平加速度峰值較小,但持續(xù)時間較長,管體的最大位移為0.02 m??傮w而言,不同地震作用下管體的位移的幅值和頻率均有差別,但主要影響因素還是地震的峰值加速度和持續(xù)時間。

        3.2 地震方向和大小

        在結(jié)構(gòu)抗震分析中,一般將地震分為水平地震作用和豎向地震作用來分析。此處通過海洋結(jié)構(gòu)抗震規(guī)范計算出豎向地震下懸浮隧道的響應(yīng)與水平地震下的響應(yīng)相比較。圖5(a)給出了不同方向地震作用下管體跨中位移時程曲線,由圖中可知水平地震的最大位移為0.30 m 方向為正,而豎直地震的最大位移為0.23 m 方向為負(fù),可以看出豎向地震產(chǎn)生的管體位移約為水平地震的76.7%。

        圖5(b)給出了不同地震方向下錨索的跨中位移時程曲線。由圖中可以看出無論是水平還是豎向地震作用下2#錨索的位移都最大。水平地震和豎向地震對應(yīng)的錨索最大位移分別為1.23 m 和0.34 m。可見水平地震會使錨索產(chǎn)生較大的位移響應(yīng),其最大位移約為豎向作用下的3.6倍。比較管體和錨索的位移,可以看出地震對于錨索的影響大于管體,這是因為錨索固定在海床上直接受到地震的作用,而管體則是通過錨索間接受到地震作用。圖5(c)給出了不同地震方向下地震峰值加速度(PGA)與隧道管體最大位移之間的關(guān)系曲線,由于采用三角級數(shù)法生成的地震譜具有一定的隨機(jī)性,系統(tǒng)響應(yīng)從而也有了隨機(jī)性使得PGA 與最大位移之間沒有形成明顯的函數(shù)關(guān)系,但對數(shù)據(jù)進(jìn)行擬合分析后發(fā)現(xiàn),PGA 與最大位移之間近似呈線性關(guān)系。

        圖4 不同地震波對懸浮隧道響應(yīng)的影響Fig.4 Impact of different seismic waves on responses of SFT

        3.3 波浪作用

        為了分析波浪作用對地震荷載作用下懸浮隧道系統(tǒng)響應(yīng)的影響,本小節(jié)分別計算了僅地震、僅波浪、地震波浪作用的三種情況。

        圖5 不同地震方向及大小對懸浮隧道系統(tǒng)的影響Fig.5 Impact of different seismic directions and PGA on responses of SFT system

        圖6(a)中給出了三種荷載情況下的管體彎矩時程曲線,可以看出地震作用和波浪作用下隧道響應(yīng)具有不同的作用特點,僅地震、僅波浪、地震波浪共同作用下最大彎矩分別為12.04×108N·m、0.18×108N·m 和15.58×108N·m,在地震作用上考慮波浪后管體彎矩增大了29.4%。圖6(b)給出了三種荷載情況下的2#錨索的張力時程曲線,僅地震、僅波浪、地震波浪作用下錨索張力分別為183.27×105N、1.65×105N、261.09×105N,在地震作用上考慮波浪后錨索張力增大了42.5%,由此可以看出對于波浪荷載對懸浮隧道系統(tǒng)的響應(yīng)具有一定促進(jìn)作用。此處僅波浪作用下錨索的位移小于僅地震作用。然而在后續(xù)計算中發(fā)現(xiàn)不同波浪參數(shù)對系統(tǒng)響應(yīng)影響較大,并非恒小于地震作用。

        進(jìn)一步的,圖6(c)~圖6(d)對波浪荷載的參數(shù)進(jìn)行了討論,選擇了不同波高、波長和周期進(jìn)行計算,為了避免地震荷載的干擾,在計算中僅考慮了波浪荷載。圖6(c)給出了波高、波長與管體最大位移的關(guān)系曲線,由圖中可以看出管體最大位移與波長、波高直接近似呈線性關(guān)系。進(jìn)一步的,圖6(d)給出了不同波浪周期下管體最大位移和波浪頻率,由圖6(d)可知隨波浪周期的增大管體最大位移出現(xiàn)明顯的峰值,最大位移達(dá)0.74 m。隨波浪周期的增大波浪頻率呈反比例下降,在波浪周期較小時(小于10 s),波浪周期與管體1階頻率相近能夠產(chǎn)生共振。

        圖6 不同波浪參數(shù)對懸浮隧道響應(yīng)的影響Fig.6 Impact of different wave parameters on responses of SFT

        4 結(jié)論

        通過Stokes三階波浪方程及三角級數(shù)法計算出波浪和地震荷載,提出了一種計算波浪地震耦合作用下懸浮隧道系統(tǒng)響應(yīng)的簡化分析模型,對荷載參數(shù)和系統(tǒng)參數(shù)進(jìn)行討論得到以下結(jié)論:

        (1)本文提出的波浪地震耦合作用下懸浮隧道簡化分析模型與錨索振動模型在計算錨索位移時具有較好的一致性,但后者不能考慮管體和錨索的耦合振動效應(yīng)。

        (2)地震的方向?qū)τ谙到y(tǒng)響應(yīng)具有顯著影響,在相同峰值加速度情況下水平地震會使系統(tǒng)產(chǎn)生更大位移響應(yīng),同期錨索的位移遠(yuǎn)大于管體。隨地震荷載峰值加速度的增加系統(tǒng)的最大位移響應(yīng)約呈線性增加趨勢。

        (3)在地震荷載的基礎(chǔ)上,考慮波浪荷載后系統(tǒng)的響應(yīng)出現(xiàn)一定增長。隨波高和波長的增大,系統(tǒng)響應(yīng)整體呈線性增大趨勢。此外當(dāng)波浪周期較小時(小于10 s)與懸浮隧道管體1階頻率相近容易引發(fā)系統(tǒng)的共振。

        猜你喜歡
        管體波浪錨索
        一種外定型無支架扣模保溫鋼管
        鋼管(2023年2期)2023-08-06 05:33:48
        波浪谷和波浪巖
        趙莊煤業(yè)13102巷圍巖補(bǔ)強(qiáng)支護(hù)技術(shù)研究與應(yīng)用
        煤(2021年10期)2021-10-11 07:06:12
        一種新型內(nèi)漲式橡膠密封圈
        波浪谷隨想
        2AI2硬鋁合金管體精密鍛造成形工藝研究
        去看神奇波浪谷
        煤礦深部高地壓巷道錨索支護(hù)技術(shù)對策研究
        不同錨索參數(shù)對預(yù)應(yīng)力場影響的數(shù)值模擬研究
        一種深孔螺釘防倒管
        科技資訊(2016年9期)2016-05-14 00:55:35
        国产裸拍裸体视频在线观看| 女人高潮久久久叫人喷水| 男人添女人囗交做爰视频| 四川发廊丰满老熟妇| 亚洲一本到无码av中文字幕| 欧美黑人乱大交| 亚州无线国产2021| 亚洲熟妇无码av不卡在线播放| 亚洲中文字幕在线一区二区三区| 欧美一级视频在线| 国产精品亚洲最新地址| 一二三区亚洲av偷拍| 欧美 日韩 人妻 高清 中文| 国产无遮挡又黄又爽在线观看| 69精品丰满人妻无码视频a片| 亚洲熟妇乱子伦在线| 日韩中文字幕乱码在线| 亚洲国产精品国自产拍性色 | 风韵犹存丰满熟妇大屁股啪啪| 国产无遮挡aaa片爽爽| 欧美成人一区二区三区| 99福利网| 在线不卡中文字幕福利| 黄片免费观看视频播放| 中文字幕人成人乱码亚洲av| 蜜桃日本免费看mv免费版| 亚洲天堂99| 97久久综合区小说区图片区| 国产诱惑人的视频在线观看| 尹人香蕉久久99天天拍| 国偷自产一区二区免费视频| 亚洲熟女少妇一区二区| 日本最新一区二区三区免费看| 国产美女高潮流白浆视频| 80s国产成年女人毛片| 国产如狼似虎富婆找强壮黑人| 人妻久久999精品1024| 精品国产日产av在线| 青春草免费在线观看视频| 在线 | 一区二区三区四区| 欧美国产日产一区二区|