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

    基于等幾何裁剪分析的拓撲與形狀集成優(yōu)化

    2016-01-07 00:58:16傅曉錦,龍凱,周利明
    振動與沖擊 2015年7期
    關鍵詞:靈敏度分析拓撲優(yōu)化有限元分析

    第一作者傅曉錦男,博士,教授,1964年11月生

    基于等幾何裁剪分析的拓撲與形狀集成優(yōu)化

    傅曉錦1,龍凱2,周利明1,闕春蘭1,葉航1

    (1. 上海電機學院機械學院,上海200245; 2. 華北電力大學新能源電力系統(tǒng)國家重點實驗室,北京102206)

    摘要:提出了將設計和分析、拓撲與形狀優(yōu)化集成的思想,探索了基于等幾何裁剪分析的拓撲與形狀集成優(yōu)化設計算法,該方法統(tǒng)一了結構優(yōu)化的計算機輔助設計、計算機輔助工程分析和優(yōu)化設計的模型,基于B樣條的等幾何裁剪分析既能準確表達幾何形狀,又可以用裁剪面分析方便處理任意復雜拓撲優(yōu)化問題,由裁剪選擇標準確定合理的拓撲結構變動方向,結構變動時無需重新劃分網(wǎng)格,設計結果突破初始設計空間的限制,還可方便優(yōu)化形狀。建立了等幾何裁剪靈敏度分析的計算方法,給出了等幾何裁剪分析拓撲與形狀集成優(yōu)化算法,通過典型實例表明所用方法的正確性和有效性。

    關鍵詞:形狀優(yōu)化; 靈敏度分析; 拓撲優(yōu)化; 等幾何分析; 裁剪B樣條; 有限元分析

    基金項目:國家自然科學基金(11202078);上海市教委科研創(chuàng)新重點項目(13ZZ145);上海市自然科學基金(11ZR1413800)

    收稿日期:2014-02-24修改稿收到日期:2014-04-21

    中圖分類號:O241;TH122

    文獻標志碼:A

    DOI:10.13465/j.cnki.jvs.2015.07.026

    Abstract:An idea of integration of computer aided design, computer aided engineering analysis, topology and shape optimization design was introduced and a kind of optimization strategy for integation of topology and shape optimization design of continuum structure based on isogeometric trimmed surface analysis was explored. The proposed method offers a unified model of computer aided design, computer aided engineering and optimization design in structure optimization. The isogeometric trimmed surface analysis based on B splines has the capability of both expressing the geometry shape accurately and solving arbitrarily complex topology optimization problem. The trimmed criteria selected determine the reasonable direction of topology changes. It does not need to remesh during structural optimization process. The design results can be beyond the limitations of initial design space. It is also convenient in shape optimization. In the integration method, the calculation approach for sensitivity analysis of isogeometric trimmed surface, and the algorithm for integation of topology and shape optimization design based on isogeometric trimmed analysis were developed. A numerical example illustrates the correctness and effectiveness of the method.

    Integration of topology and shape optimization design of continuum structure based on isogeometric trimmed surface analysis

    FUXiao-jin1,LONGKa2,ZHOULi-ming1,QUEChun-lan1,YEHang1(1. School of Mechanical Engineering, Shanghai Dianji University, Shanghai 200245,China 2.State Key Laboratory for Alternate Electrical Power System with Renewable Energy Sources, North China Electric Power University, Beijing 102206,China)

    Key words:shape optimization; sensitivity analysis; topology optimization; isogeometric analysis; trimmed B-spline; finite element analysis

    計算機輔助結構優(yōu)化設計是一項多學科相交叉的新技術。它的研究領域涉及到力學、機械、結構工程、數(shù)學和信息技術,具有巨大的發(fā)展?jié)摿1-3]。我國程耿東院士和郭旭在布局優(yōu)化研究中提出求解奇異最優(yōu)解的松弛算法,是具有重要意義的貢獻[4]。Michael等[5]提出一種的水平集方法(Level set method),優(yōu)化問題中結構的邊界用嵌入到高維尺度函數(shù)中的水平集模型來表示,模型在描述復雜結構的拓撲及邊界變化方面具有較好的靈活性。目前,多數(shù)研究者先計算機輔助設計(Computer Aided Design,CAD),之后,進行計算機輔助工程(Computer Aided Engineering,CAE)分析,幾何和分析模型會有出入,CAD/CAE模型難以融合[6-7],優(yōu)化方法的最優(yōu)解與初始產(chǎn)生空洞的數(shù)量與位置有很大關系[8-9]。其中存在二個明顯的問題:第一,設計空間與初始網(wǎng)格有關,設計空間為固定網(wǎng)格,設計結果與設計空間相關,若限制設計空間,有可能得不到解。另一方面,若網(wǎng)格劃分過細,會使計算量大大增加。第二,將計算結果返回到CAD系統(tǒng)的后處理工作量也大,這是由于常規(guī)優(yōu)化方法設計和分析的模型不是CAD系統(tǒng)中樣條數(shù)據(jù),需要大量的后處理工作。第一個缺點可采用無網(wǎng)格解決??朔诙€缺點需要統(tǒng)一設計和分析的模型[10-11]。Hughes等[12]提出了一種等幾何分析的新型有限元方法,該方法采用Non-Uniform Rational B-Splines (NURBS)作為有限元離散的插值函數(shù),將控制頂點作為網(wǎng)格節(jié)點,使CAD和CAE之間通信既方便容易,又有效[13-17]。在Hughes等[18-25]幾何分析的啟發(fā)下,很多學者在結構優(yōu)化領域結合等幾何分析進行了廣泛的研究,優(yōu)化技術迅速發(fā)展,創(chuàng)造了不少成果,具有較好的參考價值。但以上文獻的拓撲優(yōu)化方法多采用細胞或單元,然而,細胞表示的拓撲是不規(guī)則或模糊邊界輪廓,所以,試圖研究一種能確定的光滑的材料邊界的方法,最近,水平集拓撲優(yōu)化已被廣泛的研究[26],水平集拓撲優(yōu)化中,由雅可比方程確定進化水平集函數(shù)表示材料邊界,設計靈敏度分析計算形狀速度,一般,水平集拓撲優(yōu)化有一個嚴重的不足,在優(yōu)化過程中,新的內(nèi)部空洞不能自動產(chǎn)生,所以,初始的水平集會包含一些用戶定義的內(nèi)部空洞,這樣,該拓撲優(yōu)化方法的問題是會依靠初始空洞的數(shù)量和位置,在優(yōu)化過程中,合理的策略應該是基于拓撲靈敏度和應變能量密度引入新的空洞。本文在前期集成優(yōu)化研究的基礎上[27-28],探索等幾何裁剪靈敏度分析方法,利用裁剪面分析任意復雜拓撲結構,無需將其分成多塊,從而降低由于分塊引起的巨大計算量,有限元模型與幾何實體重合,統(tǒng)一了CAD、CAE和優(yōu)化設計的模型。

    1幾何造型

    非均勻有理B樣條(Non-Uniform Rational B-Splines,NURBS)為曲線曲面幾何造型主要方法[29-31]。

    1.1節(jié)點矢量與B樣條

    節(jié)點是參數(shù)軸上的分割點,一個非遞減的參數(shù)坐標的集合定義為節(jié)點矢量,記為

    Ξ=[ξ0,ξ1,…ξn+p]

    式中,ξi∈R是第i個節(jié)點,p是B樣條基函數(shù)的次數(shù),n是基函數(shù)的數(shù)量。

    定義B樣條基函數(shù)

    (1)

    (2)

    用n個p次B樣條基函數(shù)Ni,p(ξ)(i=1,2,…,n)和n個控制頂點Pi(i=1,2,…,n)確定p次B樣條曲線:

    (3)

    1.2NURBS描述

    (4)

    NURBS曲線可表示為

    (5)

    由節(jié)點矢量ξ=[ξ0,ξ1,…ξn+p]和節(jié)點矢量η=[η0,η1,…ηm+q],n×m個控制頂點Pi,j(i=1,2,…,n;j=1,2,…,m),構成的B樣條曲面:

    (6)

    由節(jié)點矢量ξ=[ξ0,ξ1,…ξn+p]、節(jié)點矢量η=[η0,η1,…ηm+q]和ζ=[ζ0,ζ1,…ζm+q],n×m×l個控制頂點Pi,j,k(i=1,2,…,n;j=1,2,…,m;k=1,2,…,l)構成的B樣條體:

    S(ξ,η,ζ)=

    (7)

    1.3基于NURBS的等幾何分析

    建立準確的幾何模型是等幾何分析的關鍵,CAD(計算機輔助設計)和CAE(計算機輔助工程分析)是同一個模型,即等幾何,用有限元方法進行分析。用B樣條基函數(shù)代替多項式基函數(shù),控制頂點視同有限元節(jié)點,可得:

    (8)

    式中,u為如位移、應力和應變等控制變量的物理場,d可以為控制頂點的位移、應力和應變等。

    2裁剪面分析

    常規(guī)的等幾何分析時,分析模型是分解成幾個B樣條小塊,所以,若拓撲結構復雜,分析較為困難,若借助計算機圖形學的裁剪操作,即使復雜模型也容易實現(xiàn)。裁剪單元按照單元的頂點數(shù)量分為三類,如圖1所示。復雜的裁剪單元與這三類裁剪單元不匹配,那么繼續(xù)四分細化該單元,直到滿足這三類的其中一類為止。裁剪單元繼續(xù)分為三角形和曲邊三角形,曲形三角形有一邊為B樣條裁剪曲線。

    圖1 三類裁剪單元 Fig.1 Three kind of trimmed elements

    因參數(shù)空間比物理空間簡單,故單元分解在未裁剪的參數(shù)空間操作。單元剛度矩陣為

    (9)

    B=MΓG

    (10)

    其中,Jξ為由參數(shù)空間向物理空間轉(zhuǎn)換的變換矩陣,表示為

    (11)

    J=JTkJTn=

    (12)

    如圖2所示,在有限元中,線性變換Tkt能用線性三角形形狀函數(shù)表示。

    圖2 三角形單元變換 Fig.2 Transformations in triangular elements

    為方便積分,采用線性變換Tkt

    Tkt:{u,v}→{ξ,η}

    (13)

    如圖3所示,正積分三角形單元的二頂點(ξ1,η1)和(ξ2,η2)與裁剪矩形單元的頂點一致,另一點(ξ3,η3)位于裁剪曲線上,可求裁剪曲線交點及其參數(shù)。假定,點(ξ3,η3)中,ξ3為已知,η3為未知。

    (1)由已知參數(shù)用下式可求λ0

    (2)由下式求得未知參數(shù)

    圖3 在三角形單元中求裁剪曲線交點及其參數(shù) Fig.3 Finding intersection and its corresponding parameter of trimming curve in triangular elements

    正三角形單元的雅可比矩陣為

    (14)

    圖4為由高斯積分單元轉(zhuǎn)化為物理空間曲邊單元,物理空間單元用ΩPH表示;相應參數(shù)空間單元用ΩPM表示。

    圖4 曲邊三角形積分單元變換 Fig.4 Transformations in curved triangular integration elements

    (1)由已知參數(shù)用下式可求λ1

    (2)由下式求得未知參數(shù)

    假定,點(ξ3,η3),ξ3為已知,η3為未知。

    (1)由已知參數(shù)用下式可求λ2

    (2)由下式求得未知參數(shù)

    圖5 在曲邊三角形單元中求裁剪曲線交點及其參數(shù) Fig.5 Finding intersection and its corresponding parameter of trimming curve in triangular elements

    在圖4中,線性變換Ti,Tj和Tkt為

    (16)

    (17)

    在αβ坐標系φ(λ)=[φα(λ),φβ(λ)]T裁剪曲線能用逆矩陣R-1變換

    φ(λ)=R-1·C(λ)=A-1(C(λ)-Z)

    (18)

    曲邊積分三角形單元的雅可比矩陣為

    J=JTiJTjJTktJTn=

    (19)

    式中,

    (20)

    圖6 積分單元分解和單元的積分點 Fig.6 Integration elements decomposition and integration points in the elements

    圖6顯示積分單元分解和積分點。單元按圖5規(guī)則分解,在曲邊三角形積分單元中的積分點均布在邊界,曲邊三角形積分單元的剛度矩陣由式(19)和式(9)計算可得。

    3裁剪靈敏度分析

    本文取邊界控制頂點的坐標為設計變量,建立優(yōu)化模型

    MinΨ(X,u(X))

    s.t.hj(X,u(X))=0j=1,…,p1

    (21)

    3.1靈敏度計算

    目標函數(shù)方程

    (22)

    式中,F(xiàn)為線性方程的力矢量,u為位移矢量。

    目標函數(shù)對設計變量的靈敏度為

    (23)

    對式(9)的單元剛度矩陣求導

    (24)

    對B的導數(shù)

    (25)

    對于矩形單元:

    對于正三角形單元:

    對于曲邊三角形單元:

    (26)

    用式(25)和式(26)代入式(24)求得單元剛度敏度,由此可得總剛度敏度

    (27)

    式中,ne為單元數(shù)量。

    形狀優(yōu)化時希望結構的最大應力最小,取最大Von-Mises應力為研究對象

    (28)

    Von-Mises應力對控制頂點坐標的導數(shù)為

    (29)

    計算極值點,求得最大應力

    (30)

    式中,d為在u全局內(nèi)最大應力點的控制頂點位移,A為由全局位移矢量轉(zhuǎn)為局部位移矢量的變換矩陣,A在該處值為1,其它地方為0。局部處靈敏度為

    (31)

    式中,v為下列伴隨方程的解。

    (32)

    4優(yōu)化方法

    本部分介紹優(yōu)化方法和策略,確定設計變量,定義選擇標準,提出相應算法。

    4.1設計變量

    在二維優(yōu)化問題,首先,未裁剪表面的邊界表面控制頂點坐標為設計變量,因為內(nèi)部表面控制頂點不影響邊界形狀,所以,內(nèi)部表面控制頂點不作為設計變量。一旦空洞生成,裁剪曲線的控制頂點坐標也作為設計變量了??刂祈旤cP的法向量n定義

    (33)

    式中,(ξP,ηP)是參數(shù)空間最接近裁剪點的控制頂點的坐標,t為切向量,k為單位向量。目標函數(shù)設計靈敏度與設計表面控制頂點法向移動有關,目標函數(shù)設計靈敏度為

    (34)

    4.2產(chǎn)生空洞的標準

    對于式(21)的優(yōu)化模型,定義選擇標準為目標函數(shù)對約束的靈敏度

    (35)

    體積作為約束條件時,SC通常為正。當空洞產(chǎn)生時,體積就會變小。如果SC為負,沒有改變約束,不可能再使目標函數(shù)性能提高。

    SC有二種,一種是設計模型內(nèi)產(chǎn)生空洞的選擇標準SC洞,另一種是設計模型內(nèi)外邊界形狀改變的選擇標準SC形。目標函數(shù)的形狀靈敏度和設計邊界的控制頂點約束都要計算相應的SC形。比較二個選擇標準SC形和SC洞能確定設計模型內(nèi)部空洞是否產(chǎn)生。

    若設計模型中,約束不變,設計控制頂點的變動而引起邊界攝動,內(nèi)部空洞產(chǎn)生,則

    (36)

    (37)

    由式(36)和式(37)可得

    (38)

    圖7 在邊界變動時表面設計控制頂點和選擇標準 Fig.7 Surface design control points and their selestion criteria on boundary movement

    4.3算法

    (1)定義模型

    (2)裁剪面分析

    (3)靈敏度分析

    (5)更新設計變量

    (6)內(nèi)部空洞是否合并?若合并返回(2)

    (7)裁剪面分析

    (8)收斂否?若不收斂,返回(3)

    (9)收斂得優(yōu)化結果

    5應用實例

    實例1

    如圖8所示為多載荷梁的結構優(yōu)化計算模型,結構尺寸為200 mm×100 mm,彈性模量E=2.01×105MPa,泊松比ν=0.3,定義整個區(qū)域為設計域。在下端分別施加100 N、600 N、100 N三個集中力載荷,并在左下角約束了6個自由度,右下角約束了y方向的移動自由度。

    圖8 優(yōu)化模型 Fig.8 Optimization model

    ●表示設計面的控制頂點,■表示初始面的控制頂點 圖9 物理空間初始面的控制頂點和設計面的控制頂點 Fig.9 Initial control pointsand design surface control points in physical area

    如圖9所示優(yōu)化模型為二次B樣條面,分為20×10個單元。在物理空間,以體積作為約束,以小方塊表示初始表面的控制頂點,以小圓點表示設計面邊界(頂端邊界)的控制頂點,首先以設計面邊界的控制頂點為設計變量,為了避免網(wǎng)格扭曲變形,二末端控制頂點允許沿垂直方向移動,其余控制頂點允許沿法向移動,其余邊界無設計面控制頂點,但可由裁剪曲線按需要改變邊界。按照上述優(yōu)化方法進行優(yōu)化,進行靈敏度分析,設計面邊界的控制頂點(即設計變量)就會按計算的靈敏度方向移動,其余表面控制頂點會相應變動,盡管在小的區(qū)域控制頂點會較密集,但自適應使其不會相交。通過改變表面控制頂點坐標進行優(yōu)化,直到由拓撲導數(shù)確定內(nèi)部產(chǎn)生空洞位置為止,得如圖10所示物理空間中表面的優(yōu)化結果。

    圖10 物理空間中表面的優(yōu)化結果 Fig.10 Optimal result in physical area

    一旦空洞生成,裁剪曲線控制頂點的坐標也成為設計變量。下一步重點分析裁剪曲線控制頂點的靈敏度。

    裁剪曲線方程為

    (39)

    式中,λ為曲線參數(shù),ng為裁剪曲線控制頂點數(shù)量,裁剪曲線第k個控制頂點坐標為Pk=[pk,qk]。

    圖11 交點靈敏度計算時裁剪曲線 控制頂點變動及其相應曲線參數(shù) Fig.11 Movement of trimming curve control point for sensitivity compute of intersections and their corresponding curve parameters

    圖11說明由于裁剪曲線控制頂點攝動,二個交點及其相應參數(shù)變動情況。攝動前,二交點分別為(ξ2b,η2b)和(ξ3b,η3b),其相應曲線參數(shù)分別為λ2b和λ3b。裁剪曲線控制頂點垂直攝動后,二交點分別變?yōu)?ξ2a,η2a)和(ξ3a,η3a),其相應曲線參數(shù)也分別變?yōu)棣?a和λ3a。交點及其參數(shù)對曲線控制頂點攝動的靈敏度分別為

    (40)

    應變位移矩陣B對曲線控制頂點的靈敏度用Γ和G的靈敏度表示為

    (41)

    其中,Γ的靈敏度表示為

    (42)

    式中,

    (43)

    G對曲線控制頂點的靈敏度表示為

    (44)

    第k個基函數(shù)對曲線控制頂點的一階導數(shù)用參數(shù)(ξ,η)靈敏度和基函數(shù)的二階導數(shù)表示

    (45)

    (46)

    由式(14)可得雅可比靈敏度

    (47)

    由式(17)計算曲線三角形積分單元的(ξ,η)靈敏度

    (48)

    由式(40)計算交點的靈敏度

    (49)

    由式(16)計算式(48)中參數(shù)α和β的靈敏度

    (50)

    由式(18)可得

    (51)

    其中

    (52)

    由式(39)可得裁剪曲線靈敏度

    (53)

    由式(15)可得曲線參數(shù)λ靈敏度

    (54)

    (55)

    (56)

    (57)

    (58)

    (59)

    對裁剪曲線求偏導

    (60)

    它對曲線控制頂點的靈敏度為

    (61)

    由式(56)到式(61)可計算曲邊三角形積分單元的式(55)雅可比行列式靈敏度。

    圖12 在參數(shù)空間中的模型Fig.12Modelintheparametricarea圖13 在物理空間中裁剪曲線控制頂點Fig.13Controlpointsinthephysicalarea圖14 在物理空間中選定裁剪曲線的詳細信息Fig.14Detailedinformationofselectcurveinthephysicalarea

    圖15 物理空間中未裁剪表面的優(yōu)化結果Fig.15Optimaluntrimmedresultinthephysicalarea圖16 物理空間表面控制頂點和逐步細化積分單元Fig.16Controlpointsofdecomposedintegrationelementsinthephysicalarea圖17 裁剪面分析應力云圖(Von-Mises應力場)Fig.17Stresscontoursoftrimmingcurve

    圖18 裁剪優(yōu)化結果Fig.18Optimaltrimmedresult圖19 形狀優(yōu)化結果Fig.19Shapeoptimalresult圖20 形狀優(yōu)化應力云圖(Von-Mises應力場)Fig.20Stresscontoursofshapeoptimization

    圖16所示為在物理空間中的模型,顯示適應表面控制頂點的敏度,圖17所示為裁剪面分析應力云圖(Von-Mises應力場),圖18所示為裁剪優(yōu)化結果。

    然后,再進行形狀優(yōu)化靈敏度分析,形狀優(yōu)化時希望結構的最大應力最小,取最大Von-Mises應力為研究對象,求Von-Mises應力對控制頂點坐標的導數(shù),計算極值點,得優(yōu)化形狀,圖19 所示為形狀優(yōu)化結果,圖20形狀優(yōu)化應力云圖。

    實例2

    如圖21所示為單載荷梁的結構優(yōu)化計算模型,結構尺寸為200 mm×100 mm,彈性模量E=2.01×105MPa,泊松比ν=0.3,定義整個區(qū)域為設計域。在右邊緣下端施加一個100 N的載荷,并約束了左邊緣6個自由度。

    如圖22所示優(yōu)化模型為二次B樣條面,分為20×10個單元。在物理空間,以體積作為約束,以小方塊表示初始表面的控制頂點,以小圓點表示設計面邊界(頂端邊界)的控制頂點,首先以設計面邊界的控制頂點為設計變量,為了避免網(wǎng)格扭曲變形,二末端控制頂點允許沿垂直方向移動,其余控制頂點允許沿法向移動,其余邊界無設計面控制頂點,但可由裁剪曲線按需要改變邊界。按照上述優(yōu)化方法進行優(yōu)化,進行靈敏度分析,設計面邊界的控制頂點(即設計變量)就會按計算的靈敏度方向移動,其余表面控制頂點會相應變動,盡管在小的區(qū)域控制頂點會較密集,但自適應使其不會相交。通過改變表面控制頂點坐標進行優(yōu)化,直到由拓撲導數(shù)確定內(nèi)部產(chǎn)生空洞位置為止,空洞生成后,裁剪曲線控制頂點的坐標也成為設計變量,更新設計變量,外部邊界和內(nèi)部邊界同時優(yōu)化,物理空間和參數(shù)空間的拓撲結構同時相應改變。

    圖21 優(yōu)化模型Fig.21Optimizationmodel圖22 物理空間初始面的控制頂點和設計面的控制頂點Fig.22Initialcontrolpointsanddesignsurfacecontrolpointsinphysicalarea圖23 在參數(shù)空間中的模型Fig.23Modelintheparametricarea

    圖24 在物理空間中裁剪曲線控制頂點 Fig.24 Control points in the physical area

    圖25 在物理空間中選定裁剪曲線的詳細信息Fig.25Detailedinformationofselectcurveinthephysicalarea圖26 物理空間表面控制頂點和逐步細化積分單元Fig.26Controlpointsofdecomposedintegrationelementsinthephysicalarea圖27 裁剪面分析應力云圖(Von-Mises應力場)Fig.27Stresscontoursoftrimmingcurve

    圖26所示為在物理空間中的模型,顯示適應表面控制頂點的敏度,圖27所示為裁剪面分析應力云圖(Von-Mises 應力場),圖28所示為裁剪優(yōu)化結果。

    然后,再進行形狀優(yōu)化靈敏度分析,形狀優(yōu)化時希望結構的最大應力最小,取最大Von-Mises應力為研究對象,求Von-Mises應力對控制頂點坐標的導數(shù),計算極值點,得優(yōu)化形狀,圖29 所示為形狀優(yōu)化結果,圖30形狀優(yōu)化應力云圖。

    圖28 裁剪優(yōu)化結果Fig.28Optimaltrimmedresult圖29 形狀優(yōu)化結果Fig.29Shapeoptimalresult圖30 形狀優(yōu)化應力云圖(Von-Mises應力場)Fig.30Stresscontoursofshapeoptimization

    6結論

    (1)本文探討了將設計和分析、拓撲與形狀優(yōu)化雙集成的基于等幾何裁剪分析的連續(xù)體結構拓撲與形狀集成優(yōu)化設計方法,該方法將在CAD系統(tǒng)表達幾何模型控制頂點和基函數(shù)的樣條信息用于CAE,CAD、CAE和優(yōu)化設計是同一個模型,實現(xiàn)其無縫集成。

    (2)用B樣條和裁剪曲線表示設計模型的內(nèi)外邊界,用裁剪面分析能方便處理復雜樣條面,改變拓撲結構,設計結果不限定在初始的設計空間內(nèi)。將未裁剪表面控制頂點坐標為設計變量,一旦空洞生成,裁剪曲線控制頂點坐標也作為設計變量,無需額外的設計變量參數(shù)。為了避免未裁剪面的幾個網(wǎng)格發(fā)生畸變,將除了角點外的設計面控制頂點沿法線方向移動,確定產(chǎn)生合理裁剪曲線的選擇標準使拓撲結構往更合理方向變動,結構變動時無需重新劃分網(wǎng)格。

    (3)建立了等幾何裁剪曲線和表面的設計控制頂點靈敏度計算方法,優(yōu)化過程中,內(nèi)部的拓撲結構是靈活可變的,提出了內(nèi)部空洞產(chǎn)生和合并算法,常規(guī)的優(yōu)化方法,設計空間相關是一個不足,提出的樣條優(yōu)化很好解決了這個問題,在拓撲優(yōu)化中,對常規(guī)網(wǎng)格相關載荷也較難處理,用提出的方法容易處理,在數(shù)字分析和優(yōu)化中,轉(zhuǎn)到CAD模型時,樣條信息是一樣的,消除后處理的工作量,同時優(yōu)化拓撲結構和形狀,典型實例表明所用方法的合理性和有效性。該方法容易擴展設計空間,方便分析處理,高度集成優(yōu)化設計,高效可行,具有工程應用價值。

    參考文獻

    [1]Huang X, Xie Y M. Evolutionary Topology Optimization of Continuum Structures: Methods and Application[M]. London: Springer, 2010.

    [2]Bendsoe M P, Sigmund O. Topology optimization: theory, methods, and applications[M]. New York:Springer-Verlag, 2003.

    [3]Qian X, Sigmund O.Isogeometric shape optimization of photonic crystals via Coons patches[J]. Computer Methods in Applied Mechanics and Engineering, 2011,200(25-28):2237-2255.

    [4]程耿東,許林. 基于可靠度的結構優(yōu)化的序列近似規(guī)劃算法[J].計算力學學報,2006,23(6):641-646.

    CHENG Geng-dong, XU Lin.Sequential approximate programming approach to reliability based structural optimization[J].Chinese Journal of Computational Mechanics 2006,23(6):641-646.

    [5]Wang M Y, Wang X, Guo D. A level set method for structural topology optimization[J]. Computer Methods in Applied Mechanics and Engineering, 2003, 192(1-2): 227-246.

    [6]Uhm T K, Kim K S, Seo Y D, et al. A locally refinable T-spline finite element method for CAD/CAE integration[J]. Structural Engineering and Mechanics, 2008,30(2):225-245.

    [7]Wall W A, Frenzel M A, Cyron C. Isogeometric structural shape optimization[J]. Computer Methods in Applied Mechanics and Engineering, 2008,197(33-40): 2976-2988.

    [8]Dedè L, Borden M J, Hughes T J R.Isogeometric analysis for topology optimization with a phase field model[J].Archives of Computational Methods in Engineering,2012,19(3):427-465.

    [9]Hassani B, Khanzadi M, Tavakkoli S M.An isogeometrical approach to structural topology optimization by optimality criteria[J]. Structural and Multidisciplinary Optimization, 2012,45(2):223-233.

    [10]Lipton S, Evans J A, Bazilevs Y, et al. Robustness of isogeometric structural discretizations under severe mesh distortion[J]. Computer Methods in Applied Mechanics and Engineering 2010, 199(5-8): 357-373.

    [11]Hassani B, Tavakkoli S M, Moghadam N Z.Application of isogeometric analysis in structural shape optimization[J].Scientia Iranica.2011,18(4):846-852.

    [12]Hughes T J R, Cottrell J A, Bazilevs Y.Isogemetric analysis:CAD, finite element, NURBS,exact geometry and mesh refinement[J]. Computer Methods in Applied Mechanics and Engineering,2005,194(39-41):4135-4195.

    [13]Cottrell J A, Reali A, Bazilevs Y, et al. Isogeometric analysis of structural vibrations[J]. Computer Methods in Applied Mechanics and Engineering, 2006, 195(41-43):5257-5296.

    [14]Auricchio F, Da Veiga L B, Hughes T J R, et al.Isogeometric collocation methods[J]. Mathematical Models and Methods in Applied Sciences, 2010, 20(11):2075-2107.

    [15]Bazilevs Y, Da Veiga L B, Cottrell J,et al. Isogeometric analysis: approximation, stability and error estimates for h-refined meshes[J]. Mathematical Models and Methods in Applied Sciences, 2006, 16(7):1031-1090.

    [16]Hughes T J R, Reali A, Sangalli G. Efficient quadrature for NURBS based isogeometric analysis[J]. Computer Methods in Applied Mechanics and Engineering, 2010,199(5-8):301-313.

    [17]Hughes T J R, The finite element method: linear static and dynamic finite element analysis[M]. Dover Publications, Mineola, New York:2000.

    [18]Roh H Y, Cho M, Integration of geometric design and mechanical analysis using B-spline functions on surface[J]. International Journal for Numerical Methods in Engineering. 2005, 62(14):1927-1949.

    [19]Li K, Qian X. Isogeometric analysis and shape optimization via boundary integral[J].Computer Aided Design,2011,43(11):1427-1437.

    [20]Hesch C, Betsch P.Isogeometric analysis and domain decomposition methods[J].Computer Methods in Applied Mechanics and Engineering,2012,213-216:104-112.

    [21]Qian X, Full analytical sensitivities in NURBS based isogeometric shape optimization[J]. Computer Methods in Applied Mechanics and Engineering.2010, 199(29-32):2059-2071.

    [22]Seo Y D, Kim H J, Youn S K. Shape optimization and its extension to topological design based on isogeometric analysis[J]. International Journal of Solids and Structures. 2010, 47(11-12):1618-1640.

    [23]Cea J, Garreau S, Guillaume P,et al. The shape and topological optimizations connection[J]. Computer Methods in Applied Mechanics and Engineering.2000,188(4):713-726.

    [24]Manh N D, Evgrafov A, Gersborg A R,et al.Isogeometric shape optimization of vibrating membranes[J]. Computer Methods in Applied Mechanics and Engineering, 2011,200(13-16):1343-1353.

    [25]Navarrina F, Gómez H,París J,et al.Isogeometric shape sensitivity analysis[J].WIT Transactions on the Built Environment,2012,125:119-130.

    [26]Allaire G, Jouve F, Toader A M. Structural optimization using sensitivity analysis and a level set method[J]. Journal of Computational Physics. 2004,194(1):363-393.

    [27]傅曉錦.連續(xù)體結構綜合優(yōu)化設計[J].機械工程學報,2012,48(1): 128-134.

    FU Xiao-jin. Integration optimization design of continuum structure[J].Chinese Journal of Mechanical Engineering, 2012,48(1): 128-134.

    [28]傅曉錦,盧炎麟.拓撲和形狀集成優(yōu)化設計[J].計算機集成制造系統(tǒng),2008,14(6): 1041-1048.

    FU Xiao-jin,LU Yan-lin. Integration of topology and shape optimization[J].Computer Integrated Manufacturing Systems, 2008,14(6): 1041-1048.

    [29]施法中. 計算機輔助幾何設計與非均勻有理B樣條[M].北京:高等教育出版社,2001.8.

    [30]Piegl L, Tiller W, The NURBS Book (2nd Edition)[M], Springer-Verlag Berlin Heidelberg, Germany: 1997.

    [31]朱心雄. 自由曲線曲面造型[M].北京:科學出版社,2000.

    [32]Cho S, Ha S H. Isogeometric shape design optimization: exact geometry and enhanced sensitivity[J]. Structural and Multidisciplinary Optimization. 2009, 38(1):53-70.

    [33]Choi K K, Kim N H. Structural sensitivity analysis and optimization 2:nonlinear systems and applications[M]. New York: Springer, 2005.

    猜你喜歡
    靈敏度分析拓撲優(yōu)化有限元分析
    基于負載均衡的電力通信流量分析和預測方法研究
    基于DOE的多軸轉(zhuǎn)向系統(tǒng)靈敏度分析
    配電變壓器對電壓偏差影響量化分析及對策
    自錨式懸索橋鋼箱梁頂推施工階段結構分析
    隨機振動載荷下發(fā)射裝置尾罩疲勞壽命分析
    航空兵器(2016年4期)2016-11-28 21:54:01
    電推衛(wèi)星方圓過渡結構輕量化設計
    有限元分析帶溝槽平封頭的應力集中
    飛機起落架支撐桿強度有限元分析
    科技視界(2016年18期)2016-11-03 22:31:14
    基于OptiStruct的襟翼搖臂的拓撲優(yōu)化
    科技視界(2016年22期)2016-10-18 15:37:54
    基于Morris法分析的液壓參數(shù)對互聯(lián)懸架的影響
    国产亚洲精品综合一区在线观看| 卡戴珊不雅视频在线播放| 在现免费观看毛片| 亚洲无线观看免费| 长腿黑丝高跟| 欧美国产日韩亚洲一区| 色综合亚洲欧美另类图片| 国产色婷婷99| 午夜激情欧美在线| 亚洲精品色激情综合| 一进一出抽搐gif免费好疼| 极品教师在线视频| 91久久精品国产一区二区三区| 亚洲电影在线观看av| 亚洲av成人精品一区久久| 免费人成视频x8x8入口观看| 日韩欧美三级三区| 久久久久免费精品人妻一区二区| 成人二区视频| 欧美精品国产亚洲| 男人舔奶头视频| 国产成人a区在线观看| 女人十人毛片免费观看3o分钟| 免费人成在线观看视频色| 听说在线观看完整版免费高清| 可以在线观看毛片的网站| 亚洲国产精品成人久久小说 | 亚洲av熟女| 99久久久亚洲精品蜜臀av| 亚洲欧美日韩高清在线视频| 久久久色成人| 综合色丁香网| 久久久久久九九精品二区国产| 国产在线精品亚洲第一网站| 精品免费久久久久久久清纯| 成人三级黄色视频| 亚洲av五月六月丁香网| 久久草成人影院| 嫩草影院精品99| 看非洲黑人一级黄片| 精品国内亚洲2022精品成人| 成熟少妇高潮喷水视频| 国产高清有码在线观看视频| 色综合亚洲欧美另类图片| 日本免费一区二区三区高清不卡| 搞女人的毛片| aaaaa片日本免费| 最近视频中文字幕2019在线8| avwww免费| av.在线天堂| 一本精品99久久精品77| 久久久精品欧美日韩精品| 日韩精品青青久久久久久| 天堂av国产一区二区熟女人妻| 欧美3d第一页| 最新在线观看一区二区三区| 又爽又黄无遮挡网站| 91精品国产九色| 精品一区二区三区av网在线观看| 精品人妻一区二区三区麻豆 | .国产精品久久| 亚洲经典国产精华液单| 三级男女做爰猛烈吃奶摸视频| 天天一区二区日本电影三级| 国产精品人妻久久久久久| 亚洲人成网站在线播放欧美日韩| 日产精品乱码卡一卡2卡三| 亚洲av中文av极速乱| 日韩欧美一区二区三区在线观看| 国产精品日韩av在线免费观看| videossex国产| 在线免费观看不下载黄p国产| 在线免费观看不下载黄p国产| 熟女人妻精品中文字幕| 久久韩国三级中文字幕| 深夜精品福利| 国产伦精品一区二区三区视频9| 搡女人真爽免费视频火全软件 | 国内少妇人妻偷人精品xxx网站| 亚洲中文日韩欧美视频| 国产色爽女视频免费观看| 日本三级黄在线观看| 亚洲一区高清亚洲精品| 国产欧美日韩一区二区精品| 久久九九热精品免费| 亚洲婷婷狠狠爱综合网| 久久精品综合一区二区三区| 久久久成人免费电影| 亚洲真实伦在线观看| 精品一区二区三区人妻视频| 亚洲无线观看免费| 最近最新中文字幕大全电影3| 亚洲久久久久久中文字幕| 国产单亲对白刺激| 人妻少妇偷人精品九色| 一区二区三区免费毛片| 色噜噜av男人的天堂激情| 亚洲最大成人中文| 大又大粗又爽又黄少妇毛片口| 午夜精品国产一区二区电影 | 国内精品美女久久久久久| 日韩大尺度精品在线看网址| 一级毛片我不卡| 日韩欧美三级三区| 婷婷色综合大香蕉| 久久久国产成人精品二区| 99久久精品热视频| 搡女人真爽免费视频火全软件 | 亚洲,欧美,日韩| 亚洲高清免费不卡视频| 少妇人妻一区二区三区视频| 亚洲中文日韩欧美视频| 国产色爽女视频免费观看| 欧美一级a爱片免费观看看| 精品日产1卡2卡| 国产精品美女特级片免费视频播放器| 久久九九热精品免费| 国产成人a∨麻豆精品| 日本免费a在线| 99久久精品热视频| 又粗又爽又猛毛片免费看| 日韩欧美在线乱码| 国产精品野战在线观看| 九九爱精品视频在线观看| 老司机福利观看| 精品日产1卡2卡| 一本一本综合久久| 午夜激情福利司机影院| 久久久色成人| 亚洲精品影视一区二区三区av| 欧美日韩乱码在线| 国产高清激情床上av| 精品国产三级普通话版| 久久热精品热| 黄色日韩在线| 成人亚洲欧美一区二区av| 国产精品国产高清国产av| 国产老妇女一区| 婷婷亚洲欧美| 国产亚洲精品av在线| 真人做人爱边吃奶动态| 亚洲av不卡在线观看| 亚洲一区高清亚洲精品| 在线观看午夜福利视频| 六月丁香七月| 国产精品久久久久久av不卡| 久久精品影院6| 成人亚洲欧美一区二区av| 老司机影院成人| 一个人免费在线观看电影| 一级毛片久久久久久久久女| 乱人视频在线观看| 少妇裸体淫交视频免费看高清| 亚洲国产欧美人成| 久久久久久久午夜电影| 亚洲精品国产成人久久av| 在线免费观看的www视频| 免费电影在线观看免费观看| 亚洲精品日韩在线中文字幕 | 国产精品伦人一区二区| 国产中年淑女户外野战色| 欧美成人精品欧美一级黄| 三级毛片av免费| 熟女人妻精品中文字幕| 简卡轻食公司| 99久久无色码亚洲精品果冻| 国产成人福利小说| 亚洲乱码一区二区免费版| 国产女主播在线喷水免费视频网站 | 亚洲人成网站在线播| 国产片特级美女逼逼视频| 国模一区二区三区四区视频| 欧美xxxx黑人xx丫x性爽| 亚洲va在线va天堂va国产| 99热只有精品国产| 日本与韩国留学比较| 久久热精品热| 亚洲第一区二区三区不卡| 白带黄色成豆腐渣| 波野结衣二区三区在线| 精品不卡国产一区二区三区| 人人妻人人澡人人爽人人夜夜 | 国产伦精品一区二区三区视频9| 亚洲欧美精品自产自拍| 国产精品不卡视频一区二区| 黄色欧美视频在线观看| 国产在线精品亚洲第一网站| 看片在线看免费视频| 十八禁国产超污无遮挡网站| 老师上课跳d突然被开到最大视频| 国产精品永久免费网站| 午夜激情欧美在线| 99热全是精品| 国产高清有码在线观看视频| 六月丁香七月| 在线免费观看不下载黄p国产| 狂野欧美白嫩少妇大欣赏| 欧美成人免费av一区二区三区| 又爽又黄a免费视频| 最新在线观看一区二区三区| 国产精品久久久久久久久免| 三级毛片av免费| 男插女下体视频免费在线播放| 卡戴珊不雅视频在线播放| 一级毛片久久久久久久久女| 91在线精品国自产拍蜜月| 欧美丝袜亚洲另类| 成人综合一区亚洲| 中国美白少妇内射xxxbb| 不卡一级毛片| 搡老岳熟女国产| 国产一区二区三区在线臀色熟女| 成人性生交大片免费视频hd| 国产精品伦人一区二区| 精品久久久久久久久久久久久| 可以在线观看毛片的网站| 国产白丝娇喘喷水9色精品| 精品熟女少妇av免费看| 国产精品福利在线免费观看| 亚洲av免费在线观看| 成人性生交大片免费视频hd| 老师上课跳d突然被开到最大视频| 亚洲,欧美,日韩| 天天躁夜夜躁狠狠久久av| 热99在线观看视频| 九九热线精品视视频播放| 日日干狠狠操夜夜爽| 春色校园在线视频观看| 干丝袜人妻中文字幕| 国产aⅴ精品一区二区三区波| 国产成人一区二区在线| 亚洲激情五月婷婷啪啪| 久久这里只有精品中国| 在线观看美女被高潮喷水网站| 午夜激情欧美在线| 精品久久久久久久久久久久久| 国产午夜精品论理片| 日韩一本色道免费dvd| 国产精品免费一区二区三区在线| 男女之事视频高清在线观看| 老司机影院成人| 国产精品伦人一区二区| 看片在线看免费视频| 中文字幕免费在线视频6| 久久久久精品国产欧美久久久| 国产精品不卡视频一区二区| 欧美三级亚洲精品| 国产精品无大码| 长腿黑丝高跟| 18禁裸乳无遮挡免费网站照片| 午夜福利成人在线免费观看| av专区在线播放| 欧美激情久久久久久爽电影| 亚洲成人精品中文字幕电影| www.色视频.com| 深爱激情五月婷婷| 18禁在线无遮挡免费观看视频 | 免费在线观看成人毛片| 亚洲欧美精品综合久久99| 国产精品乱码一区二三区的特点| 国产高清激情床上av| videossex国产| 欧美色欧美亚洲另类二区| 亚洲天堂国产精品一区在线| 欧美激情国产日韩精品一区| 国产精品久久久久久亚洲av鲁大| 五月伊人婷婷丁香| 欧美zozozo另类| 欧美极品一区二区三区四区| 亚洲无线观看免费| 免费看a级黄色片| 少妇丰满av| 哪里可以看免费的av片| 日韩av在线大香蕉| 精品人妻熟女av久视频| 成人欧美大片| 熟女电影av网| 久久这里只有精品中国| 亚洲av第一区精品v没综合| 日韩 亚洲 欧美在线| 亚洲在线观看片| 成人高潮视频无遮挡免费网站| 精品人妻偷拍中文字幕| 精品国产三级普通话版| 中国国产av一级| 国产黄色视频一区二区在线观看 | 特级一级黄色大片| 日韩欧美精品v在线| a级毛片a级免费在线| 美女 人体艺术 gogo| 天堂网av新在线| 在线国产一区二区在线| 99视频精品全部免费 在线| 国产一区二区三区av在线 | 伦理电影大哥的女人| 三级经典国产精品| 国产成人精品久久久久久| avwww免费| 亚洲精品国产成人久久av| 高清毛片免费观看视频网站| 精品久久久久久久久久久久久| 亚洲最大成人中文| 午夜福利在线观看免费完整高清在 | 亚洲成人久久爱视频| 中出人妻视频一区二区| 午夜免费男女啪啪视频观看 | 97在线视频观看| 午夜老司机福利剧场| 欧洲精品卡2卡3卡4卡5卡区| 亚洲性久久影院| 国产午夜福利久久久久久| 日韩欧美一区二区三区在线观看| 久久精品综合一区二区三区| 国产成人影院久久av| 亚洲欧美精品综合久久99| 尾随美女入室| 国产亚洲精品av在线| 亚洲色图av天堂| 高清午夜精品一区二区三区 | 精品国内亚洲2022精品成人| 成年版毛片免费区| .国产精品久久| 蜜臀久久99精品久久宅男| 美女黄网站色视频| 国产一区二区亚洲精品在线观看| 狠狠狠狠99中文字幕| 国产精品嫩草影院av在线观看| 夜夜夜夜夜久久久久| 男女那种视频在线观看| 人人妻人人澡人人爽人人夜夜 | 又黄又爽又免费观看的视频| 国产精品人妻久久久久久| 在线播放无遮挡| 搡老熟女国产l中国老女人| 永久网站在线| 午夜a级毛片| 长腿黑丝高跟| 美女cb高潮喷水在线观看| 午夜福利在线观看吧| 最近的中文字幕免费完整| 99热这里只有是精品50| 亚洲性久久影院| 在线观看av片永久免费下载| 18禁黄网站禁片免费观看直播| 亚洲av美国av| 美女xxoo啪啪120秒动态图| 精品人妻视频免费看| av中文乱码字幕在线| 露出奶头的视频| 国产片特级美女逼逼视频| 国产高清三级在线| 精品久久久久久久人妻蜜臀av| 中文亚洲av片在线观看爽| 五月伊人婷婷丁香| av在线天堂中文字幕| 午夜精品在线福利| 色尼玛亚洲综合影院| 亚洲人成网站在线观看播放| 人人妻人人澡人人爽人人夜夜 | 久久精品人妻少妇| 91av网一区二区| a级毛片免费高清观看在线播放| 毛片女人毛片| 日本三级黄在线观看| 特级一级黄色大片| 日韩精品青青久久久久久| 免费不卡的大黄色大毛片视频在线观看 | 久久精品国产自在天天线| 91麻豆精品激情在线观看国产| 国产一区二区在线av高清观看| 国产精品久久久久久av不卡| 欧美成人精品欧美一级黄| 18+在线观看网站| 免费在线观看影片大全网站| 日韩三级伦理在线观看| 亚洲精品久久国产高清桃花| 亚洲图色成人| 白带黄色成豆腐渣| 国产成人91sexporn| 日韩成人av中文字幕在线观看 | 久久久久国产精品人妻aⅴ院| 日本色播在线视频| 精品免费久久久久久久清纯| 亚洲国产精品久久男人天堂| 国语自产精品视频在线第100页| 黑人高潮一二区| 国产精品久久久久久久久免| 午夜免费男女啪啪视频观看 | 免费大片18禁| 国产欧美日韩一区二区精品| 欧美日韩国产亚洲二区| 日产精品乱码卡一卡2卡三| 一个人免费在线观看电影| 特大巨黑吊av在线直播| 网址你懂的国产日韩在线| eeuss影院久久| 午夜老司机福利剧场| 观看免费一级毛片| 1000部很黄的大片| 国产精品女同一区二区软件| 色av中文字幕| 亚洲中文字幕一区二区三区有码在线看| 久久精品国产清高在天天线| 中国美女看黄片| 在线国产一区二区在线| 美女大奶头视频| 搡老妇女老女人老熟妇| 精品午夜福利视频在线观看一区| 免费观看在线日韩| АⅤ资源中文在线天堂| 乱系列少妇在线播放| 内射极品少妇av片p| 99热只有精品国产| 午夜精品在线福利| 欧美日韩国产亚洲二区| 菩萨蛮人人尽说江南好唐韦庄 | 精品欧美国产一区二区三| 联通29元200g的流量卡| av中文乱码字幕在线| 性插视频无遮挡在线免费观看| 久久久久久久久久久丰满| 欧美日韩一区二区视频在线观看视频在线 | 在线观看av片永久免费下载| 偷拍熟女少妇极品色| 久久草成人影院| 婷婷精品国产亚洲av| 丰满人妻一区二区三区视频av| 乱系列少妇在线播放| 国产成人freesex在线 | 免费在线观看成人毛片| 波野结衣二区三区在线| 亚洲国产欧美人成| 欧美另类亚洲清纯唯美| 不卡一级毛片| 亚洲av第一区精品v没综合| 老师上课跳d突然被开到最大视频| 久久欧美精品欧美久久欧美| 一边摸一边抽搐一进一小说| 免费人成在线观看视频色| 国产一区二区三区av在线 | 亚洲性夜色夜夜综合| 内地一区二区视频在线| ponron亚洲| 丝袜美腿在线中文| 欧美日韩国产亚洲二区| 亚洲美女视频黄频| 国产单亲对白刺激| 亚洲激情五月婷婷啪啪| 亚洲四区av| 久久久久国产网址| 国产精品国产三级国产av玫瑰| 最新在线观看一区二区三区| 99热这里只有是精品在线观看| 听说在线观看完整版免费高清| 极品教师在线视频| 日韩精品青青久久久久久| 国产亚洲av嫩草精品影院| 国产成人freesex在线 | 午夜精品国产一区二区电影 | 在线观看美女被高潮喷水网站| 国内精品久久久久精免费| 精品久久久久久久久亚洲| 国产欧美日韩一区二区精品| 免费av不卡在线播放| 噜噜噜噜噜久久久久久91| 国产精品永久免费网站| 99国产极品粉嫩在线观看| 成人无遮挡网站| 在线播放无遮挡| 免费不卡的大黄色大毛片视频在线观看 | 黄片wwwwww| 亚洲国产色片| 国产精品一区二区三区四区免费观看 | 久久精品综合一区二区三区| 校园春色视频在线观看| 12—13女人毛片做爰片一| 久久韩国三级中文字幕| 日韩成人伦理影院| 国产欧美日韩精品亚洲av| 97超级碰碰碰精品色视频在线观看| 全区人妻精品视频| 亚洲av一区综合| 国内久久婷婷六月综合欲色啪| 欧美日韩在线观看h| 乱系列少妇在线播放| 波野结衣二区三区在线| 69人妻影院| 午夜日韩欧美国产| 亚洲综合色惰| 成人鲁丝片一二三区免费| 22中文网久久字幕| 国产精品久久久久久av不卡| 欧美日韩国产亚洲二区| 变态另类丝袜制服| 色尼玛亚洲综合影院| 精品久久久久久久久亚洲| 日本黄色片子视频| 身体一侧抽搐| 国产成人a∨麻豆精品| 日韩成人伦理影院| 长腿黑丝高跟| 国产精品永久免费网站| 2021天堂中文幕一二区在线观| 久久午夜亚洲精品久久| 天天一区二区日本电影三级| 噜噜噜噜噜久久久久久91| 人妻丰满熟妇av一区二区三区| 亚洲美女视频黄频| 干丝袜人妻中文字幕| 日韩人妻高清精品专区| 免费观看的影片在线观看| 欧美最新免费一区二区三区| 乱人视频在线观看| 精品久久久久久久人妻蜜臀av| 在线免费观看的www视频| 一卡2卡三卡四卡精品乱码亚洲| 日日啪夜夜撸| 日韩亚洲欧美综合| 日本成人三级电影网站| 国产黄色小视频在线观看| 国产在线男女| 精品午夜福利在线看| 国产大屁股一区二区在线视频| 国产高清激情床上av| 如何舔出高潮| 18禁黄网站禁片免费观看直播| 成年免费大片在线观看| 国产三级中文精品| 亚洲国产色片| 精品熟女少妇av免费看| 六月丁香七月| 久久久久国产网址| 深爱激情五月婷婷| 久久久欧美国产精品| 久久久久久国产a免费观看| 国产成人a区在线观看| 久久久久久伊人网av| 蜜桃久久精品国产亚洲av| 亚洲精品日韩av片在线观看| 久久精品夜色国产| 精品一区二区三区视频在线观看免费| 卡戴珊不雅视频在线播放| 久久久久国产网址| а√天堂www在线а√下载| 午夜免费男女啪啪视频观看 | 久久久国产成人精品二区| 国产精品久久久久久亚洲av鲁大| 乱系列少妇在线播放| 午夜福利在线观看免费完整高清在 | 日韩欧美 国产精品| 亚洲欧美清纯卡通| 你懂的网址亚洲精品在线观看 | 亚洲av二区三区四区| 男女之事视频高清在线观看| 欧美日本视频| 亚洲熟妇熟女久久| 日韩大尺度精品在线看网址| 亚洲自偷自拍三级| 久久久a久久爽久久v久久| 久久中文看片网| 观看免费一级毛片| 国产精品国产三级国产av玫瑰| 日韩欧美国产在线观看| 久久久成人免费电影| 校园春色视频在线观看| 精品久久久久久久久亚洲| 极品教师在线视频| 亚洲成人中文字幕在线播放| 又爽又黄无遮挡网站| 亚洲自拍偷在线| 99久国产av精品| 国产一区二区在线观看日韩| 欧美另类亚洲清纯唯美| 国产一区二区在线观看日韩| 成人亚洲精品av一区二区| 日本色播在线视频| 久久草成人影院| 少妇人妻一区二区三区视频| 成人鲁丝片一二三区免费| 搡老岳熟女国产| 女的被弄到高潮叫床怎么办| 久久精品久久久久久噜噜老黄 | 成人av一区二区三区在线看| 精品欧美国产一区二区三| 精品一区二区三区人妻视频| 1024手机看黄色片| 亚洲成人久久性| 国产美女午夜福利| av天堂在线播放| 中出人妻视频一区二区| 一区福利在线观看| 在线播放国产精品三级| 亚洲综合色惰| 国产色爽女视频免费观看| 日韩欧美三级三区| 99热这里只有是精品50| 国产高清视频在线观看网站| 别揉我奶头~嗯~啊~动态视频| 身体一侧抽搐| 乱系列少妇在线播放| 亚洲av中文av极速乱| 最新在线观看一区二区三区| 一级a爱片免费观看的视频| 99久久成人亚洲精品观看| 亚洲av免费在线观看| 美女被艹到高潮喷水动态| 亚洲国产欧美人成| 欧美日韩国产亚洲二区| 两个人视频免费观看高清| 黄色日韩在线| 久久久色成人| 日日啪夜夜撸| 免费看日本二区| 丰满的人妻完整版| 成人综合一区亚洲|