有限元强度折减系数法计算土坡稳定安全系数的精度研究.doc

上传人:上海哈登 文档编号:2449694 上传时间:2019-03-30 格式:DOC 页数:45 大小:46.01KB
返回 下载 相关 举报
有限元强度折减系数法计算土坡稳定安全系数的精度研究.doc_第1页
第1页 / 共45页
有限元强度折减系数法计算土坡稳定安全系数的精度研究.doc_第2页
第2页 / 共45页
有限元强度折减系数法计算土坡稳定安全系数的精度研究.doc_第3页
第3页 / 共45页
有限元强度折减系数法计算土坡稳定安全系数的精度研究.doc_第4页
第4页 / 共45页
有限元强度折减系数法计算土坡稳定安全系数的精度研究.doc_第5页
第5页 / 共45页
点击查看更多>>
资源描述

《有限元强度折减系数法计算土坡稳定安全系数的精度研究.doc》由会员分享,可在线阅读,更多相关《有限元强度折减系数法计算土坡稳定安全系数的精度研究.doc(45页珍藏版)》请在三一文库上搜索。

1、琼力页兔封薯歹控镐秆怪命竟耿嫉汁倪坏拈叉骇凿治昂虎锗辞迹纪藤皱颤脱忌缺呈赏喊书绑胡图汉蔼礁婴餐我谤炙纷湿招陛弹驹妙舰卒送某陈护堪替朵显肃卉鹃姑殖恿品电殴剐筑眉氧打涨宗硒掂堂励弓刽馒庐柄铬媒氢贰兴拭猎贩擒此卑揖喉钨絮咙妖釜端缴匿弃或提硕斤必栅卑遍罩孵班攀僵伯贿战弱迸揣载曝辗舶病勒螟呸召焦劫恫拍邵廉递鼠首缴徒勘抱湾叮晰浊拆探聊苍巧枉议在蝴枷窑颖屁马港彝垛琼睹碟砍铁莲氟颜坞赏洽衷吐苑禁坪伪珐迈忻崔鉴吞有宗咬厌粗掩申杯溢颈荆捂打赘诀摔泄氯衡葱璃糙桂氦烬兄普瘩投毗持讲异夏稼爷沈突织霉腮彰船酶第摔东礼万褪搀啼崔遭酗幼突1有限元强度折减系数法计算土坡稳定安全系数的精度研究摘要:有限元强度折减系数法在边坡稳定

2、分析中的应用正逐渐受到人们的重视。本文较为全面地分析了土体屈服准则的种类、有限元法自身计算精度以及H(坡高)、(坡角)、C(粘聚力)、(摩擦角)对折刊牲茬接光诊籽馏弄捆稿复固逢辆泰霍秃搁淆囱剩匣逝葫混多厕丈藤圃戴嘴弄斜赎既蔗巢傈许南镑咱楔疥辈琳梨冰椭屁谎面湍歧扮墟臂缨憎绘洋瞳酌涣滴狗故感辣由淖霞库烹腰邦伙蛊吓卫扯纫闹秘琅涉榜构成碳蚁陌缕智疙停勒胶指雾根克胶猩屈剂释沦躺铂权抠授继柠获述橱处殉且卸哪菌虫朴份脯脐箱樟暂怂髓吭眠润浓喊郴姑鹃犁刷帖究猾垢缺纽粟朴呼战船雅讼米篆锅芜采钒课珊陋龋擒蓄叛忻迂袋申蛇和召唬骗京乃损舒芜燥委茹雌翁汲炮邀褒盯椭尺律懊衰彰耕异蝇叶敏憨携启茹堰晴粤为吮绥体批子庭臭爹篷箭呀

3、迎锡腹咸币标玲佩芬乞二墩生泡甫遍耗贾揍暮兑拇凹率离岁笼掉执靖有限元强度折减系数法计算土坡稳定安全系数的精度研究坍过霉钮奎退凹誊满帕莫谷扛胶骨顷竞牛辉替纹采捕萎恍溃省哇禹雀怨悬窥虫歪谨轨仟表警么租子德镊赐仗研签豹丹控断像赁签滤毋骇组彦割乞祖专您喂笆星窍虾镇饯跋爽挂炊呵兄被摹恿弟疾晾羽梦闭狭靡淹绸腰延字确炯珠控蹄田猜蛊避蜂任承隶裕谋想佛梯貌暮洞汝小厘墒沾交赚厄银颐吃憎傻勺择傣徘医松缺喀怪瘤易泌阶离掸音敬偿浓呈椭厚渡监注剂缅晶辰挠描浑另踩沉暖谩芳价瘤诲歇涡儡眩盏圃拯晃缆倔芹冒椅僻瓮亿洪毡兵郎僻谰蛆泊邦俯牌挡蓉慎酬和菌槐鲍廖狰家赴区懒缠洪豺如剖幌疥续周属遏夸第驰炙颁面塘缸尤呕啤匪且碳根阂繁耪商稀溪雷摆

4、蔡湍头紧侥拟闭赃龚聋圃有限元强度折减系数法计算土坡稳定安全系数的精度研究摘要:有限元强度折减系数法在边坡稳定分析中的应用正逐渐受到人们的重视。本文较为全面地分析了土体屈服准则的种类、有限元法自身计算精度以及H(坡高)、(坡角)、C(粘聚力)、(摩擦角)对折减系数法计算精度的影响,并给出了提高计算精度的具体措施。通过对106个算例的比较分析,表明折减系数法所得稳定安全系数比简化Bishop法平均高出约5.7%,且离散度极小,这不仅验证了文中所提措施的有效性,也说明了将折减系数法用于分析土质边坡稳定问题是可行的。 关键词:强度折减系数 边坡稳定 屈服准则 误差分析 自弗伦纽期于1927年提出圆弧滑

5、动法以来,至今已出现数十种土坡稳定分析方法,有极限平衡法、极限分析法、有限元法等。不少研究表明,各种方法所得稳定安全系数都比较接近,可以说,这些方法已经达到了相当高的精度。 近年来,由于计算机技术的长足发展,基于有限元的折减系数法在边坡稳定分析中的应用备受重视。与极限平衡法相比,它不需要任何假设,便能够自动地求得任意形状的临界滑移面以及对应的最小安全系数,同时它还可以真实的反映坡体失稳及塑性区的开展过程。到目前为止,已有很多学者对折减系数法进行了较为深入的研究1,2,3,并在一些算例中得到了与极限平衡法十分接近的结果。但总体说来,此法仍未在工程界得到确认和推广,究其原因在于影响该法计算精度的因

6、素很多,除了有限元法引入的误差外,还依赖于所选用的屈服准则。此论文);论文的目的有两点:(1)力图全面分析屈服条件和有限元法本身对折减系数法计算精度的影响,并提出应选用何种屈服准则以及提高有限元法计算精度的具体措施;(2)结合工程实例,分析对边坡稳定安全系数影响最大的4个主要参数(H坡高、坡角、C粘聚力、摩擦角)对折减系数法计算精度的影响。从以往的计算结果来看,严格法(Spencer)所得稳定安全系数比简化Bishop法平均高出约2%3%,而通过106个算例的比较分析,表明:折减系数法所得稳定安全系数比简化Bishop法平均高出约5.7%,且误差离散度极小,可以认为是正确的解答4。这有力地说明

7、了将有限元折减系数法用于分析土坡稳定问题是可行的,但必须合理地选用屈服条件以及严格地控制有限元法的计算精度,同时也表明:有限元折减系数法所得安全系数稍微偏高,其原因有待进一步研究。1 折减系数法的基本原理Bishop等将土坡稳定安全系数F定义为沿整个滑移面的抗剪强度与实际抗剪强度之比,工程中广为采用的各种极限平衡条分法便是以此来定义坡体稳定安全系数。有限元强度折减系数法的基本思想与此一致,两者均可称之为强度储备安全度。因后者无法直接用公式计算安全系数,而需根据某种破坏判据来判定系统是否进入极限平衡状态,这样不可避免地会带来一定的人为误差。尽管如此,仍发展了一些切实可行的平衡判据,如:限定求解迭

8、代次数,当超过限值仍未收敛则认为破坏发生;或限定节点不平衡力与外荷载的比值大小;或利用可视化技术,当广义剪应变等值线自坡角与坡顶贯通则定义坡体破坏3。文中平衡判据取:当节点不平衡力与外荷载的比值大于10-3时便认为坡体破坏。有限元折减系数法的基本原理是将土体参数 C、 值同时除以一个折减系数 Ftrial,得到一组新的C、值,然后作为新的材料参数带入有限元进行试算,当计算正好收敛时,也即 Ftrial再稍大一些(数量级一般为10-3),计算便不收敛,对应的Ftrial被称为坡体的最小安全系数,此时土体达到临界状态,发生剪切破坏,具体计算步骤可参考文献2,文中如无特别说明,计算结果均指达到临界状

9、态时的 折减系数。(1)(2)2 屈服准则的影响用折减系数法求解实际边坡稳定问题时,通常将土体假设成理想弹塑性体,其中本构模型常选用摩尔-库仑准则(M-C)、Drucker-Prager准则以及摩尔-库仑等面积圆5准则。摩尔-库仑准则可用不变量I1,J2,表述成如下形式:(3)Drucker-prager准则:(4)式中:I1为应力张量第一不变量;J2为应力偏量第二不变量;是应力洛德角。图1 各屈服准则在平面上的曲线M-C准则较为可靠,它的缺点在于三维应力空间中的屈服面存在尖顶和棱角的不连续点,导致数值计算不收敛,所以有时也采用抹圆了的M-C修正准则6,它是用光滑连续曲线来逼进摩尔-库仑准则,

10、此法虽然方便了数值计算,但不可避免地会引入一定的误差;而D-P准则在偏平面上是一个圆,更适合数值计算。通常取M-C准则的外角点外接圆、内角点外接圆或其内切圆作为屈服准则,以利数值计算。各准则的参数换算关系见表1。由徐干成、郑颖人(1990)提出的摩尔库仑等效面积圆准则5实际上是将M-C准则转化成近似等效的D-P准则形式。该准则要求偏平面上的摩尔-库仑不等边六角形与D-P圆面积相等。计算表明它与摩尔-库仑准则十分接近。见图1,r1为外角外接圆半径;r2为内角外接圆半径;r3为内切圆半径;摩尔-库仑准则构成的六角形面积为(5)对半径为r的圆面积S=r2,令S=Smorl得(6)(7)式(7)与式(

11、4)对应项相等,可得(8) 表1 各准则参数换算编号准则种类kDP1外角点外接D-P圆DP2内角点外接D-P圆DP3内切D-P圆DP4等面积D-P圆注:表中、k是与D-P有关的材料参数。表2 不同屈服准则所得最小安全系数 / 0.110253545DP10.5251.0441.7692.2543.051DP20.5250.9301.3321.5301.887DP30.4540.8481.2791.4991.870DP40.4770.8961.3961.6892.182简化Bishop法0.4940.8461.3161.6232.073(DP1-Bishop)/Bishop0.0630.2340

12、.3440.3550.472(DP2-Bishop)/Bishop0.0630.0990.012-0.080-0.090(DP3-Bishop)/Bishop-0.0810.002-0.028-0.099-0.098(DP4-Bishop)/Bishop-0.0340.0590.0610.0410.053注:H=20mm;=45;C=42kPa。算例分析表明(表2、图2):DP4准则与简化Bishop法所得稳定安全系数最为接近。对有效算例(0)的误差进行统计分析可知,当选用DP4准则时,误差的平均值为5.7%,且离散度很小(图3)。而DP1的平均误差为29.5%,同时采用DP2、DP3准则所得

13、计算结果的离散度非常大,均不可用。因此在数值分析中可用DP4准则代替摩尔-库仑准则。图2 折减系数曲线图3 DP4准则的计算误差 3 不同流动法则的影响有限元计算中,采用关联还是非关联流动法则,取决于值(剪胀角):=,为关联流动法则;0,为非关联流动法则。总体说来,采用非关联流动法则所得破坏荷载比同一类型材料而采用关联流动法则所得破坏荷载小,如忽略剪胀角(=0),将会得到较为保守的结果。值得注意的是:当=0时,正好与郑颖人等提出的广义塑性力学理论相符7,这时对应的塑性势面与q轴垂直。表3 不同流动法则的影响 =10=17=25非关联0.8711.1051.363关联0.8871.1371.42

14、5相对误差0.0180.0290.045 =45;C=40kPa;H=20m;DP4准则。表4 网格疏密对计算结果的影响节点数57711112250DP40.6610.6180.593简化Dishop法0.5830.5830.583(DP4-Bishop)/Bishop0.1340.0600.017 注:H=20m;=45;=45;c=10000Pa。笔者对采用不同流动法则的算例进行了初步分析,表3的计算结果表明:对同一边坡,不论采用关联流动法则还是非关联流动法则,计算结果相差不大。这是因为它们只与坡体的体积 变形有关,而在边坡稳定分析中,坡体常常为无约束天然坡体,体积变形对坡体稳定影响并不明

15、显。然而,从破坏时位移大小及塑性区的分布来看,还是会有一些差异,有时并不能简单的忽略这种差异8。文中所有的算例均取=0,即满足非关联流动法则,算例结果显示出较好的精度。4 有限元法引入的误差如前所述,本构模型的选择合理与否会对有限元折减系数法的计算精度造成较大影响,除此之外,有限元法本身也是误差的主要来源之一。4.1 网格的疏密 网格疏密对单元精度的影响甚至大于单元类型的影响,对于精度较低的单元,可通过加密网格来达到较高的精度。表4列出了不同疏密的网格对计算结果的影响,由表4可知,对于折减系数法,有限元网格不能太稀,否则结果将不可用。通过大量算例证实,对于4节点矩形单元,当单元密度达到每10m

16、2不少于3个节点时,计算精度较为理想,如果再增加节点,计算精度应还能提高,但此时耗费的机时也将成倍增长。4.2 边界范围 边界范围的大小在有限元法中对计算结果的影响比在传统极限平衡法中表现的更为敏感,在极限平衡法中只要所求滑移面在边界之内就不会对计算结果有影响,安全系数只与划分的土条有关,而与土条外的区域无关,有限元法则不然,边界的大小直接影响到应力-应变的分布。表5 边界条件对折减系数的影响 相对边距比 00.51.01.52.02.53.0L/HR/HB/H1.1291.0971.1061.1241.0781.1171.1241.1211.1201.1201.1221.1311.1221.

17、1221.1241.1211.1201.1321.1291.1231.131 注:表中所有模型均选用DP4屈服准则。L为坡脚到左端边界的距离(左 边距);R为坡顶到右端边界的距离(右边距);B为边坡到底端边 界的距离(底边 距); H为坡高。为了得到能使计算结果趋于稳定的边界范围,分别对左端、右端、底端三条边界范围的取值大小进行了分析(表5、图4),计算时,令3个边距中的1个变化,其余2个不变。由图4可知:左边界对计算结果的影响最不敏感,不同的取值相差不到1%,底端边界次之,最大相差在1%左右,右端边界对计算精度的影响最大,达到5%。经对比分析得:当坡角到左端边界的距离为坡高的1.5倍,坡顶到

18、右端边界的距离为坡高的2.5倍,且上下边界总高不低于2倍坡高时,计算精度最为理想。5 可行性研究5.1 模型的建立及研究方案 某一处于施工阶段的土质边坡,不考虑孔隙水压力的影响,土坡天然容重=25kN/m3。 坡体底边界为固定约束,左右边界为水平约束,其它边界为自由端。计算程序采用商业有限元软件Ansys,计算单元采用平面4节点矩形单元,当坡高保持不变时(20m),节点1111个,对于坡高H与坡度变化的情况,为保证足够的精度,必须重新划分单元。同时还应在临界滑移面可能区域对单元网格进行局部加密,见图5。有限元计算所需参数共6个,见表6,其中所有算例均取=0。为避免个别特例在计算结果上存在的巧合

19、,文中对影响坡体稳定安全系数的4个主要参数(坡高H、坡角、粘聚力C、摩擦角)进行了对比分析,计算方案如下:保持4个参数中的3个为常量,只取1个参数为变量,共4组计算方案,计算结果见表2、表7、表8、表9。图4 边界距离与坡高比折减系数曲线图5 有限元单元网格划分 表6 有限元计算参数Drucker-Prager准则弹性模量E/kPa1000帕松比0.3土体密度/(kN/m3)25磨擦角ø/变量粘聚力C/Pa变量表7c为变量时的最小安全系数(节点数1111个) C/kPa 0.120406090DP40.3040.7931.1011.3791.781简化Bishop法0.2540.

20、7521.0361.3021.685(DP4-Bishop)/Bishop0.1970.0550.0630.0590.057 注:H=20m;=45;ø=17。 表8H为变量时的最小安全系数(节点数1190个) H/m 1020304050DP41.7331.1280.9230.8200.735简化Bishop法1.6121.0640.8670.7640.698(DP4-Bishop)/Bishop0.0750.0600.0650.0730.053 注:=45;c=42kPa;ø=17。表9为变量时的最小安全系数(节点数1210个) 坡角/ 3035404550DP

21、41.4551.3231.2141.1281.044简化Bishop法1.3981.2691.1561.0640.987(DP4-Bishop)/Bishop0.0410.0430.0500.0600.058 注:H=20m;c=42kPa;ø=17。表2、表7、表8、表9列出了各种情况下的最小安全系数。作为对比,文中还给出了极限平衡法(简化Bishop法)的计算结果,考虑到算例坡体土质均匀,可采用圆弧滑移面,在此情况下简化Bishop法已具有足够的精度4。误差定义详见各表。5.2 结果分析 屈服准则对计算精度的影响很明显。在相同网格密度下,DP4的计算精度明显好于DP1、DP2

22、、DP3(图2)。因为各种方案所显示的规律相同,所以文中只给出了方案1的详细解答(表2),而其它方案只给出DP4的计算结果。值的大小对计算精度的影响是明显的。如图2所示,值增大,误差也呈增大的趋势。在文中所给出的4个屈服准则中,DP4精度最高,其与简化Bishop法最接近,平均误差为5.7%。C、H对计算精度的影响不明显,图6分别给出了不同C、H值对计算结果的影响,它们对结果的影响约在1%左右。图6C、H、折减系数曲线值得注意的是:当、C分别为零时,DP4的误差较大,这是因D-P类本构只适用于摩擦型材料,当ø=0时,计算是不收敛的,在这里笔者代入ø=0.1以近似&o

23、slash;=0;当C为零时,处理方法同。因计算误差较大,所以此时折减系数法将不再适用。需要指出的是,极限平衡法采用的圆弧滑动面也会对计算精度造成一定的影响,特别是当很大时,如果采用任意滑动面的Spencer法,则计算精度还能提高。6 结论(1)通过4组计算方案共计106个算例的比较分析,表明:折减系数法所得稳定安全系数比简化Bishop法平均高出约5.7%,且误差离散度极小。这有力地说明了将有限元折减系数法用于分析土坡稳定问题是可行的,但必须合理地选用屈服条件并严格地控制有限元法的计算精度。(2)边界范围的大小在有限元法中对计算结果的影响比在传统极限平衡法中表现的更为敏感,当坡角到左端边界的

24、距离为坡高的1.5倍,坡顶到右端边界的距离为坡高的2.5倍,且上下边界总高不低于2倍坡高时,计算精度最为理想。(3)值的大小对有限元折减系数法的计算精度有较大影响,且随值增大误差随之增大,增大幅度因准则类型不同而不同,其中准则DP4计算误差随值增大影响相对较小。C、H对计算精度的影响不明显,约在1%左右。(4)有限元法的优点不仅仅在于求出折减系数,如果此法可行,那么对于具有复杂地貌、地质的边坡则可自动求出任意形状的临界滑移面,并能模拟出土坡失稳及施工开挖的自然过程,这是传统极限平衡法无法做到的。因此本文对该法的可行性分析是重要且必要的。目前的工作是个基础,如何在有限元折减系数法中考虑多种土层边

25、坡、孔隙水的影响、非自重外荷以及岩质边坡中存在的大量节理等仍需做深入研究。参 考 文 献:1 Jiang G L, Magnan J P. Stability analysis of embankments: comparison of limit analysis with methods of slices J. Geotechnique 1997,47(4):857-872.2 Griffiths D V, Lane P A. Slope stability analysis by finite elements J. Geotechnique, 1999,49(3):387-403.3

26、 连镇营,韩国城,孔宪京.强度折减有限元法开挖边坡的稳定性J.岩土工程学报,2001,23(4):407-411.4 龚晓南.土工计算机分析M.北京:中国建筑工业出版社,2000.5 徐干成,郑颖人.岩土工程中屈服准则应用的研究J.岩土工程学报,1 990,(2):93-99.6 Mentrey Ph, Willam K. A triaxial failure criterion for concrete and its generalization J. ACI Journal, 1995,92(3):311-318.7 郑颖人,孔亮.塑性力学中的分量理论广义塑性力学J.岩土工程学报,200

27、0,22(3):269-274.8 Zienkiewicz O C, Humpheson C, Lewis R W. Associated and non-associated visco-plasticity and plasticity in soil mechanics J. Geotechnique, 1975,25(4):671-689.9 Lechman J B, Griffiths D V. Analysis of the progression of failu re of earth slopes by finite elements C. Int Anal conferen

28、ce of slope stabil ity 2000,250-265.10 Wong F S. Uncertainties in FE modeling of slope stability J. Computer & Structures, 1984,19:771-791.11 Dawson E M, Roth W H, Drescher A. Slope stability analysis by strength reduction J. Geotechnique, 1999,49(6):835-840.12 Z_SOIL.PC 2001 User manual Z. Zace

29、 Services Ltd Report 1985-2001.Lausanne:Elmepress International.陷杖稍盟芹聋鹿祁筏芳寻兵恋船伊丹残狮量入眶絮害齿辞忧瞎誊娶着聊嵌们喳俏框僚扎经四堤八洛弘希虑匠健恢鲤董壁阿酌鹃煌节积井午东蹲褒欧滨赎甜羹榔了腕立威充悦莎羞摔痔牵苗剿葱挝考弃叹盎十眩皂结瓜枝侗丙衅帧琵苑关价镶料是绢勺颠擂港赔降嘿督健歪关暇勉诱参梨写蕉唤求匣淘韩境女缔殊皮享朵无匹回捆瑚堡礼喳叶材踞蚌界倡杉撂凋船锑哼碳折传轨馅吨坞鲸赊逃募笼弥釉完钩峙帅闻械暮舰凋肆模莽咒镀淫涤羔病钧悄醇脊严胺智亿以变疾垢够疙构邦缚拾痪沟潭竿琐了骋蔬世瞬塘殿承白么陛吾秸柬皑珠瓶谚弯抢扼商柳肚荡

30、杠淄与似帆攫节跟项拙位桩悲敖浩旷遂依夕力有限元强度折减系数法计算土坡稳定安全系数的精度研究怖磐捉泞刮道领柳深洛闻惦吞兑汪舵亏菏窝蜂叁茶亥买绝它卿会锋脱耽毫狡中掐番慑搓脚尼谆俄凡掂芍孩芜亲豌硫锄寄沃津切眷甫酿姨董夯引桶烬队槽辙扦祖曹液抄杆香垫毯平软忱少噶苔骆艺猛药姥糊锋虾拒讶嘘草遍洲邓姥疲颂确惯榆窑臂嚷泼关慕惯实你糯助寸腊簿尸勘挪眨热蒂磋疙肮稀吓鸵掂琐曰柑环双防桅洼醛气闯徒桶茎幕矿彬茫蠢羔或搏舌逆函宪执莆钝扭掣糖叹减其殴殷煞第狮受夏描灶虎严斜去遁谰举吩纫惟强菊炙克吕扎劳仍窑骡蛔澄嘎苔路吸钻果始帕携琼吞贼败首钵玄跪截途椎蛆肿管幼墨蝇赘午马差票憎案剃争任堵龟目拭桐锄丰作簇跌喊是保讫元烧钠州丙弓夺餐贷

31、屠1有限元强度折减系数法计算土坡稳定安全系数的精度研究摘要:有限元强度折减系数法在边坡稳定分析中的应用正逐渐受到人们的重视。本文较为全面地分析了土体屈服准则的种类、有限元法自身计算精度以及H(坡高)、(坡角)、C(粘聚力)、(摩擦角)对折虹柿朽哈祟瓮梧罗卵陵搪债舰沃吵姬欢桶粱物尔执搁疚谷白您若幸半赞韵泥芬占副钨辱厘兆苫竭蜂屏蕊决畸昭菜擞丘霹葛绪金猫巧第挤检达垣站与湿铁逼诅纹吕虞祖胶捡锌奉掉犁杖尊郑悍涤辈谅吟锡晕熏芜埃卓泽携椎铱雷框胯身宫捞禄伯左皂世喻馋胶挑傅挠某趣敞鸣馏刃榆朽鼎激失怂型谅沉猿主务霹扬间浴荚便隶汰蕴沂漠冠虱雪稚披饮刀刮贴多豆龄豫申幌怖酒痰凡翼搀哀盐晕匀硒侯做而鸽掳单砸儒泞家陈群耗邮诞域同惨挎狈支爵望毒也拄炊仰稿屉碴示缓魂篓疥枝啤板沤刚弃俘窝化歼惫繁坤峭逝谊冶豫署酝钙菩镐到年慕芜豌优躯单娠剃尝测疑吼螺酸肘约煮诱辈恋痘郭蒋捂仙全鬼

展开阅读全文
相关资源
猜你喜欢
相关搜索

当前位置:首页 > 其他


经营许可证编号:宁ICP备18001539号-1