三维有限元法

三维有限元法（共11篇）

三维有限元法 篇1

0 引言

螺栓连接是工程结构中被普遍采用的一种连接方式。作为连接件,螺栓的强度关系着包括被连接件在内的整个结构的安全性与可靠性。目前国内在螺栓连接强度设计方面还没有形成统一的规范。螺栓或螺纹连接的强度分析主要采用的方法有解析法和有限元法。解析法在螺栓强度分析计算处理过程中做了大量简化,且部分参数计算繁杂,从而限制了其广泛应用[1,2,3]。近20年来,数值计算技术的发展使有限元法逐渐成为螺栓连接强度分析的主要手段。有限元法的精度与采用的模型密切相关,模型主要包括梁单元模型[4]、二维轴对称模型[5]、三维轴对称模型[6,7]以及三维螺旋模型[8,9,10]等。梁单元模型未考虑螺纹细节,规模小;二维轴对称模型考虑螺纹细节,但忽略了螺旋效应的影响,且无法应用于偏置偏载场合;三维轴对称模型能应用于偏置偏载场合,但忽略了螺纹升角的影响;三维螺旋模型精确,但规模较大,在处理大量螺栓组时计算成本过高。文献[8-10]虽采用三维螺旋模型,但都忽略了螺纹收尾和端部倒角造成的不完整螺纹对计算精度的影响。因此,针对螺栓连接的强度分析,研究一种精度较高、建模便捷以及计算成本较低的有限元分析方法具有十分重要的意义。

本文采用有限元软件ABAQUS建立螺栓连接的精确三维螺旋模型和三维轴对称简化模型,分别获得螺栓连接的应力分布特征,并对两种建模方式下的螺栓螺纹段应力分布进行了分析比较,探讨了两种建模方法的差别。

1 螺栓有限元模型

螺栓连接系统主要包括螺栓、螺母和夹持件,见图1。本文使用三种不同的有限元模型对螺栓连接系统进行分析:三维螺旋模型采用289 190个单元和300 797个节点的六面体网格划分,单元类型为C3D8IC3D6,并以少量楔形体和四面体单元进行过渡实现精确建模,并考虑螺纹收尾和端面倒角造成的不完整螺纹的影响,见图2a;轴对称模型中采用子模型方法分析普通偏置偏载场合的螺栓强度,见图2b;三维无螺纹简化模型结合解析法进行螺栓强度设计,见图2c。本文以某公司的某动力总成悬置产品为例,螺栓的螺纹结构采用细牙螺纹,螺纹规格为ISO公制标准三角形螺纹M10×1.25,螺栓各部分尺寸参照ISO公制系列六角螺栓和六角螺母的规定值。

2 螺栓计算边界条件

2.1 预紧力的计算

螺栓预紧方式主要有扭矩控制式、转角控制式和屈服控制式三种。从设计角度,由强度等级、设计载荷来确定螺栓的最大预紧应力,根据预紧方式确定预紧载荷;从校核计算角度,根据拧紧力矩或螺母转角确定螺栓的预紧力。对屈服控制式预紧,预紧力可近似取设计手册中给出的最小拉力载荷或保证载荷。对以拧紧力矩给出的预紧载荷,拧紧力矩与预紧力的关系为

则式(1)可整理为

式中,Tf为拧紧力矩;Fpre为预紧力;P为螺纹节距;dp为外螺纹有效直径,计算时可取dp=d2;ρ′为当量摩擦角;α为牙形半角;β为螺纹升角;μs为螺纹副摩擦因数;μw为螺母支撑面摩擦因数;dw为支撑面等效直径;B为螺母对边距;Di为螺栓孔内径。

对以螺母转角给出的预紧载荷,螺母转角与预紧力的关系为[1?3]

式中,φ为螺母转角;Kt为螺栓刚度;Kc为夹持件刚度。

2.2 预紧力的计算

有限元分析中常用的预紧力施加方式有温度载荷法、过盈配合法和预紧截面法。温度载荷法把预紧载荷换算成温度梯度载荷施加到螺栓光杆部分,使螺栓在热载荷作用下产生收缩获得预应力状态。但此方法需通过试算、迭代获得温度载荷,不适用于包含有其他热载荷边界条件的应用场合;同时在研究螺栓本身的应力时,很难通过控制其径向与轴向收缩获得准确的应力。过盈配合法在建立有限元建模时设定夹持件之间的过盈配合关系产生预应力,通过试算、迭代获得过盈量并合理布置过盈接触位置,从而获得准确的螺栓应力分布,常用于研究螺栓应力。预紧截面法是多数商用有限元软件提供的一种简便的预紧力加载手段,通过将预紧力直接施加在螺栓的截面上,产生大小相等、方向相反的法向作用力获得预紧力,此时螺栓的应力分布略有偏差,而夹持件的应力分布较为准确,在研究夹持件时应用较多。

本文在分析螺栓应力分布特征时采用过盈配合法,对夹持件进行分析时采用预紧截面法施加预紧力,并就两种方法获得的螺栓应力分布进行对比,见图1。称螺栓应力分析工况为工况Ⅰ,该工况中,在夹持件A与B相互接触给定过盈量(过盈量通过试算获得,试算过程是先施加一定的过盈量,输出螺栓截面上的合力,使得该合力与螺栓预紧力相等,重复此过程,通常3次试算即可获得较好的结果,本工况中取δc=0.0382mm);称夹持件分析工况为工况Ⅱ,其中,在螺栓光杆部位某横截面上施加预紧力(Fpre=8.0kN)。

3 计算结果分析

3.1 三维螺旋模型计算结果

3.1.1 螺栓纵向应力分布

螺栓材料为30CrNiMo8(1),其机械性能参数见表1。

由螺栓连接等效应力(图3)与轴向正应力分布(图4)可见,螺栓的高应力区域主要集中在螺栓头部与杆身过渡圆角、杆身与螺纹过渡部位以及承载螺纹的前端。螺栓的光杆中段和非承载螺纹中段的应力呈近似均匀分布;螺杆头部、螺纹收尾和承载螺纹附近区域的应力等值面近似呈旋转对称的抛物面,但受螺纹升角、不完整螺纹以及承载螺纹的周向接触特性等因素的影响,连接刚度沿周向分布不均,致使应力等值抛物面最高点与螺栓轴线略有偏离。这种应力分布特征从图5所示实际螺栓断口的宏观及内部裂纹形态[11,12,13]也可明显地观察到。螺纹升角是螺栓周向刚度分布不均匀的根本原因。螺纹升角越小,不完整螺纹和承载螺纹的周向接触越趋于均匀,抛物面越接近对称形态。

夹持件的应力等值面呈以接触面为分界线的两个轴对称中空球台形态(图4),可见螺栓周向刚度的不均匀分布对夹持件压应力等值面的对称性影响较小。解析法在处理夹持件刚度时所采用的简化模型便是基于这种分布形态的,目前较为常用的是中空锥体和中空圆柱体等效方法[3]。

沿螺栓的轴向正应力与轴向位移分布见图6。在螺栓杆身与非承载螺纹段,应力为常值,螺栓头部与杆身收尾段因轴向刚度变化,应力呈现渐变趋势。在承载螺纹段应力随啮合螺纹长度的增加而减小。位移曲线呈现出明显的双线性特征,分别为杆身(k1)和非承载螺纹(k2)区域,并在刚度渐变区域和承载螺纹段呈现非线性过渡。

3.1.2 螺栓横向应力分布

螺栓应力在横截面上的分布见图7。在靠近螺纹端的横截面上受螺纹升角的影响,周向刚度分布不均匀,应力等值线分布出现偏斜;在A、C截面之间横截面上的应力偏差幅度约为0.5%。采用预紧截面法加载时,预紧截面一定程度上破坏了螺栓的连续性,从而对螺栓的应力造成一定影响。为最小化这种影响,现将预紧截面放置在A、C截面中间的B处,见图8。由图8可见,通过合理放置预紧截面,可获得与过盈配合法非常接近的螺栓应力分布,保持因刚度沿周向分布不均而造成的应力分布特征。而在靠近螺栓头部的位置,应力等值线呈现同心圆形式,由于受轴向刚度的不连续性影响,因此无论是采用过盈配合加载还是预紧截面加载,螺栓头部横截面上的应力分布均呈现出类似的不均匀性。

3.1.3 螺纹应力分布

等效应力沿螺栓外轮廓线路径方向上的分布情况见图9。Ⅰ为螺栓头部与杆身过渡圆角处的应力集中点,应力值随圆角半径增大而减小;光杆段应力维持常值至接近螺纹收尾位置,因收尾螺纹牙根造成应力释放而下降;Ⅱ为收尾螺纹牙根的应力集中点,随过渡螺纹牙根的加深应力集中呈增加趋势(Ⅲ点),在沿周向渐变时,Ⅲ点应力值将随Ⅱ点的增加而减小;非承载完整螺纹段牙根最大应力大致维持常值,在接近承载螺纹和收尾螺纹端,因应力集中的影响而略微走高;最大应力发生在第一扣承载螺纹前端的牙根部位((1)点),随后呈快速下降趋势,且下降速度逐渐减缓,直至最后一扣承载螺纹前部的牙根部位((6)点);最后一扣承载螺纹后部牙根部位的应力值比前部牙根处有所增加((7)),该分布特点与文献[5]相同。

图10所示为螺纹牙根等效应力沿螺旋线路径的分布情况。ac为收尾螺纹段,应力呈逐渐上升趋势,b点峰值是周向刚度不均匀引起的应力集中;ce为非承载完整螺纹段,应力为常值且伴有周向刚度不均匀引起的波动;de段受第一扣承载螺纹的影响,应力较ce有所增加;ef为第一扣承载螺纹前部的牙根应力,e点是接触不连续造成的跳动,由于螺母的内螺纹为不完整螺牙(端部倒角),接触刚度随螺旋角增大而增大,因而应力缓慢上升;fg为后续承载螺纹前部牙根应力,呈逐渐下降趋势;g为接触终止点,应力有跳跃;gh段为最后一扣承载螺纹后部牙根应力,其应力水平比该扣螺纹前部牙根处略高,此外,因螺母的最后一扣螺纹不完整,接触刚度变化较大而引起应力波动较为明显。

螺纹根部等效应力沿轴向分布如图11所示,M、S和E分别为牙根中点、接触开始点和接触终止点。在非承载螺纹段(t1)最大应力发生在根部圆角的中心(M点);随着内外螺纹进入啮合,由于连续两扣承载螺纹的相互影响,最大等效应力位置逐渐接近接触起始点S(t2~t6),并呈逐渐减小的趋势;最后一扣承载螺纹后部牙根的应力(t7)因无后续承载螺纹的影响,其最大应力位置向上部螺纹偏移,且比前部的牙根应力略高(t6)。此外由于边接触的应力集中,应力在S点和E点的值较高。由此可见,文献[2]进行解析分析时假设螺纹剪切断裂裂纹起始位置在接触起始点是合适的,而文献[3]的处理方法略显粗糙。

3.2 三维轴对称简化模型计算结果

三维螺旋模型精度高,但建模困难且计算成本较高,在实际工程中较少使用。轴对称模型因便于参数化建模且可考虑螺纹细节,成为当前主要的简化手段。目前广泛采用的有二维轴对称模型和三维轴对称模型。前者模型规模小,常用在螺栓受同轴载荷作用的场合;后者则常用于偏置偏载场合。

如果轴对称简化后完全忽略了螺纹升角的影响,应力等值面呈完全轴对称的抛物面形态,是一个理想的应力分布状态,然而这种情况下的计算结果偏离了实际情况。因此采用三维轴对称简化模型进行螺栓的强度计算时,应考虑螺纹升角的影响,获得的承载螺纹段等效应力的分布情况见图12。而最大主应力分布表明(图13),最大拉应力发生在第一扣承载螺牙前部牙根,在承受动载荷时,较易在此处萌生疲劳裂纹。

轴对称模型的等效应力沿外轮廓线方向、牙根等效应力沿轴向分布和轴向正应力与变形分布情况分别如图14~图16所示。通过图9与图14对比分析,结果表明,轴对称模型中承载螺纹段最大等效应力下降趋势的非线性程度比精确模型中同等条件下最大等效应力下降趋势的非线性程度小2.4%。通过图11与图15对比分析,结果表明将轴对称截面放置在不完整螺纹段的中部时,轴对称模型与精确模型中承载螺纹段各峰值应力最为接近,误差不超过1.3%,所以轴对称模型中牙根应力分布能很好地保持精确模型的牙根应力分布形态。通过图6与图16对比分析,结果表明,在轴向正应力作用下,精确模型中螺纹段轴向位移线性变化斜率在0.5%左右,而轴对称模型中螺纹段轴向位移线性变化斜率在0.56%左右,所以轴对称模型的轴向力学特性受螺纹升角影响的偏差较小。

因此,轴对称简化模型的横向力学特性呈完全对称分布,虽然与精确模型有本质区别,但是能较好地反映螺栓的轴向力学特性,仅在承载螺纹段的应力水平受螺纹升角的影响有一定差异。文献[9-10]研究表明,当螺纹升角β≤4°时,螺旋效应对螺栓应力的影响较小,轴对称模型的精度基本满足工程应用要求。

3.3 三维轴对称简化模型的应用

将三维轴对称模型与子模型求解技术相结合,可有效地提高处理大量螺栓组计算的速度。

某动力总成悬置的螺栓连接见图17,在橡胶安装圆孔位置受两个方向力,力的作用平面与橡胶轴向对称面有一定偏置。全局模型中采用简化的螺栓模型,获得承载最恶劣的螺栓位置,然后对该位置螺栓采用三维轴对称子模型进行详细分析。图18中A部分为全局模型中的螺栓模型,B部分为三维回转轴对称子模型。采用节点位置映射的方法进行计算,结果发现轴对称子模型的位移边界与全局模型的连续性保持一致,子模型的计算结果精确可行,其精度能够满足工程应用要求。

4 结论和展望

(1)本文主要讨论建模方法在螺栓连接有限元分析中的应用,对粗牙和细牙螺纹均适用。

(2)采用预紧截面加载方法时,通过合理布置预紧面的位置可获得较为精确的结果。

(3)精确的三维螺旋模型的结果显示,同轴载荷作用下,由于受螺纹升角的影响,螺栓应力、应变分布呈非对称性,但这种非对称分布特征在夹持件上表现不是很明显。

(4)由于螺纹升角的存在,在螺栓的收尾过渡螺纹以及螺母的不完整螺纹段造成轴向刚度和接触刚度沿周向分布不均匀,这是螺栓纵向力学特性的主要特性。但在螺纹升角不大的情况下,采用三维轴对称简化模型,在控制计算规模的同时得到了较精确的应力结果,通过实例进行验证,可以获得与精确三维螺旋模型较为接近的螺栓纵向力学特性。

(5)三维轴对称模型虽解决了建模困难的难题,但在处理大量螺栓组时,同样面临模型规模过大的问题。另外,文中所讨论的螺栓连接应力分析是在材料的弹性范围内,因此,在今后的研究中,将需要进一步考虑材料的塑性效应和大量螺栓组统一建模等问题。

三维有限元法 篇2

动载作用下沥青路面车辙三维有限元分析

针对沥青路面纵坡路段的车辙问题,采用矩形波动荷载加载,利用三维有限元数值法对重载,慢速及高温条件下路面各结构层竖向位移和竖向应力分布规律进行分析.结果表明:荷载越大、车速越慢、温度越高,沥青路面结构的竖向位移及竖向应力越大,沥青路面纵坡路段车辙越容易形成.

作 者：梁涛 王艳荣 Liang Tao Wang Yanrong 作者单位：陕西省铜川公路管理局,陕西,铜川,727031刊 名：现代交通技术英文刊名：MODERN TRANSPORTATION TECHNOLOGY年，卷(期)：7(1)分类号：U416.217关键词：道路工程 车辙 有限元 动荷载 影响因素

三维有限元法 篇3

关键词：ansys软件；上部结构；动力分析；反应谱

中图分类号：TV312 文献标识码：A 文章编号：1006-8937（2016）18-0064-03

1 概 述

水电站厂房是能量进行转化的场所，将水能通过工程机械转化成机械能并最终转换为电能。它通过合理的工程手段，使河水平顺地引入水轮机，能量转化后引出水轮机，同时水电站厂房也为能量转化的设备提供合适的安装位置，为这些设备的安装、检修和运行提供方便的条件。

水电站厂房结构比较特殊，厂房各构件尺寸庞大，内部各结构受力条件复杂。厂房的主要组成结构包括：上游—下游挡水墩墙、钢筋混凝土蜗壳、导叶、尾水管、上部结构。

上部结构由主厂房下游的柱墙结构和副厂房上的板梁柱结构组成。上部结构很高，而宽度和厚度相对于高度方向尺寸很小，在地震效应作用下很容易发生摆动，使结构发生破坏，进而影响厂房的整体运行，造成重大的损失，因此对上部结构的抗震分析是很有必要的。

2 计算模型和计算理论

2.1 计算模型

某水电站位于四川省境内，为二等大（2）型工程，工程正常蓄水位398 m，相应库容6 330万m3，装机容量4×190 MW+1× 12 MW（生态机组），额定水头33 m。枢纽布置包括河床式电站、船闸、13孔泄洪冲砂闸、左岸副坝、左岸非溢流坝、右岸接头坝。本文通过建立一个机组段厂房三维有限元模型来了解地震作用下现有结构布置方案结构的应力、应变。计算采用的直角坐标系为：X轴为沿水流方向，顺水流方向为正，Y轴正方向为竖直方向，向上为正，Z轴垂直于水流方向，河流右岸为正方向，坐标原点高程为327 m。整体三维有限元模型，如图1所示。

2.2 计算理论

在地震作用下，结构系统的有限元方程是

2.3 计算假定

①混凝土、基础岩体为均质、弹性、各向同性的连续体，不考虑钢筋和混凝土的应力重分布。

②厂房为多个坝段，各厂房坝段之间分别设有结构缝，在计算时，因此各坝段独立承担荷载，坝段间无相互作用。

③计算时，结构中的二期混凝不承受荷载。。

2.4 荷载及荷载组合

抗震计算水位为上游正常蓄水位、下游最低尾水位，静荷载包括厂房自重、设备自重、内水压力、侧水压力、扬压力，动水压力的影响采用目前坝工界普遍采用的韦斯特加特（Westergarrd）公式进行计算：

式中：Pw为作用在坝体单位面积上的动水压力；

αh为设计地震加速度水平向代表值；

h为计算位置距水面的深度；

H为库水总深度；

ρw为水的密度。

3 模态分析

模态分析是研究结构振动特性的方法，能够确定自振频率、机型参与系数及振型等结构的振动特性。模态分析是在进行其他动力分析之前进行的，主要是由于结构的振动特性决定结构对于各种动力荷载的响应情况。厂房整体结构的前20阶自振频率，见表1。厂房上部结构的自振频率，见表2。可以看出整体的自振频率比较密集，其中九阶表现为上部结构的自振。厂房坝段第五、六阶振型图，如图2和图3所示。

4 动应力分析

本工程所处区域的地震设计烈度为7.3 °，水平向设计地震加速度代表值αh=0.13 g，竖向设计地震加速度代表值αv= 0.087 g，设计反应谱按《水工建筑物抗震设计规范》中4.3.3采用，设计反应谱最大值的代表值βmax=2.25，最小值不应小于 βmax=0.45，场地类别为I类，响应特征周期Tg=0.20 s，由此确定抗震计算所用的设计反应谱，如图4所示。

5 动应力和位移结果分析

抗震分析时其荷载按照顺水流方向、坝轴线方向和竖直向三个方向同时受地震荷载作用，其中竖向地震荷载为水平向的2/3，结构的总动力响应为顺水流方向和坝轴线方向动力响应的平方和开平方与竖向动力响应的0.5倍直接相加；最终总地震效应为反应谱计算的地震动应力和静力计算得到的静应力的叠加。动静叠加时需要对地震作用效应按系数0.35进行折减，再与静力计算结果进行叠加。反应谱分析得到的动应力是交变应力，所以在进行动应力和静应力叠加时，应分别进行正向叠加和负向叠加。由于地震作用下可能导致止水失效，故本次计算时按照扬压力系数为0.6和1.0两种工况进行计算。地震工况的最大位移表，见表3。

①分析应力结果可知，主厂房下游柱子与发电机层的相交部位X向的拉应力约为3.5 MPa，属于体型结构突变处，产生应力集中，但应力集中范围较小。主厂房下游柱子与副厂房上部板梁柱之间的联系梁以及副厂房上部板梁柱结构与下游墩墙之间的连系梁X向拉应力很大，约为5.0 MPa。副厂房上部板梁柱结构中板柱相交处X向拉应力交大，约为3.0 MPa。上部结构中的Z向梁的Z向拉应力较大，约为3.5 MPa。

②分析位移结果可知，主厂房下有柱子顶部X向最大位移为4.47～4.72 cm；主副厂房吊车梁顶部最大位移为3.70～3.93 cm，位移均比较大，但柱子底部的X向最大位移为2.55～2.74 cm，可知其位移是由厂房的整体位移造成的。经计算，吊车梁轨顶侧向位移满足《水电站厂房设计规范SL 266-2001》中表4.2.7的要求。

5 结 语

①上部结构中，主厂房下游柱子与发电机层相交处以及副厂房上部板梁柱结构中的板柱相交部位拉应力较大，可以通过适当的配筋来提高该部位的抗拉性能，有利于结构的安全。

②主厂房下游柱子之间的纵向连系梁、主厂房下游柱子与副厂房上部板梁柱结构之间的连系梁以及副厂房上部板梁柱结构与下游尾水墩墙之间的连系梁均有很大的拉应力，一定程度上减小了上部结构在X和Z向的摆动，有助于结构抗震。

③通过对上部结构的应力和位移分析，除个别部位有较大的拉应力外，其它部位的拉应力较小，应力分布符合一般规律，满足设计要求。位移也满足规范要求，结构合理。

参考文献：

[1] 侯攀，陈尧隆，邓瞻.用ANSYS对水电站厂房坝段进行抗震分析[J].西北 水力发电，2005，（1）.

[2] SL266-2001中国人民共和国行业标准.水电站厂房设计规范[S].

[3] 刘启钊.水电站[M].北京：中国水利水电出版社，1997.

[4] 党国强，李守义，鞠静春，等.河床式水电站厂房坝段动力分析[J].电网与 水力发电进展，2008，（3）.

[5] 樊锐，陈尧隆，刘武军，等.河床式水电站厂房坝段三维有限元抗震分 析[J].水资源与水工程学报，2009，（5）.

[6] 揽生瑞，杨菊生，李守义，等.河床式水电站结构分析与结构特性研究[J].西安理工大学学报，1994，（4）.

三维有限元法 篇4

颞下颌关节是人体中最为复杂的滑膜关节之一, 为负重关节, 是具有转动和滑动的双重关节, 担负着言语、咀嚼、吞咽等功能[3]。此关节由关节窝、髁突及软骨、关节盘、关节囊以及韧带构成, 其结构复杂, 形态精细并且不规则, 组织之间性质不同, 在关节运动时, 不同的结构应力不同, 这些因素使得建立高质量的颞下颌关节模型更加困难, 故对颞下颌关节的每一个细微结构的模拟显得尤其重要, 国内外专家也一直寻求高效、准确的建模方法。本文对关于颞下颌关节每一结构的建模进行回顾。

1下颌骨有限元模型的建立

目前, 有限元方法建模主要有:磨片、切片法, 三维测量法, CT图像处理法, CT、MRI图像联合处理法, DICOM数据建模法, 数字化虚拟人建模法。

1.1磨片、切片法通过实体模型磨切的方法, 扫描每一断层, 获得相应的图像信息后, 再将图像的顺序叠加, 最后建立三维模型。此方法属于破坏性建模, 误差较大;仅适用于离体的模型;颞下颌关节结构比较复杂, 现已很少适用此方法建模。Tanaka[4]在研究紧咬牙时颞下颌关节不同区域的应力分布时用到此方法。

1.2 三维测量法 用接触或非接触方法进行测量, 把数据输入计算机, 建立三维牙颌模型。高勃等[5]用此方法进行牙冠表面测量。但由于此方法不能反映模型内部情况, 故不适用于颞下颌关节的建模。

1.3 CT图像处理法 目前研究口腔生物力学较常用的方法, 对志愿者进行CT扫描, 获得原始二维数据, 将数据输入计算机, 获得二维图像, 运用图像分析软件获取图像的边界数据, 最终将获得的数据输入三维有限元分析软件中进行处理, 获得三维有限元模型。

1.3.1 髁突表面软骨层的模拟:髁突覆盖着一层纤维软骨, 在前斜面较厚, 区别于透明软骨, 它具有缓冲、分散关节的载荷, 润滑关节、减小摩擦和磨损的作用。国内大多数专家建立髁突软骨的厚度为0.2mm, 源自Pullinger[6]所做的研究数据。 胡凯等[7,8]在研究咬合时颞下颌关节应力分布时, 建立的模型中, 关节软骨的厚度为0.2mm;孔亮等[9]建立5种不同张口位颞下颌关节模型时, 模拟的关节软骨厚度为0.2mm;杨辉、周学军、胡林华等[10~12]根据关节软骨的解剖特征, 所建立的模型中其厚度为前份0.8mm, 后份1.0mm, 这比全部厚度都为0.2mm更为精确;刘展[13]在模拟颞下颌关节内软组织时, 关节功能区 (髁突前斜面和关节结节后斜面) 的软骨层最厚, 采用0.5mm;髁突和关节窝顶部最薄, 采用0.2mm。

1.3.2 下颌骨材料属性的选择:材料参数的选择很多, Nagahara[14]在用成年干颅骨研究紧咬牙时颞下颌关节位移和应力分布时, 所用的材料属性见表1。

1.3.3 模型网格化:在建模时, 适当的网格划分可以使操作简化, 网格分割越细, 模型的几何形态越接近原物, 就能得到越来越接近精确的结果[15]。在模型网格化, 划分单元时, 单元的类型选择也很重要。一般应选择实体单元, 而且尽量少用不稳定的四面体单元。六面体单元最为合适[16]。

1.4 MRI图像处理法 MRI图像可直观地显示颞下颌关节区软组织 (如关节盘和盘附着) 结构的形态和位置[17]。可以根据扫描模型获得的信息建模, 弥补了CT扫描显示关节软骨不佳的缺陷, 提高了颞下颌关节模型的几何相似性;MRI具有软组织对比度分辨力高的特点, 对骨组织形态显示较好, 因此可以此为数据建立包括髁突、关节盘、关节窝在内的TMJ三维有限元模型[18]。

1.5 DICOM数据建模法 患者完成CT扫描后, 不必生成胶片, 大量的数据信息可依照DICOM标准存盘或传输。运用相应的软件, 如Mimics等, 输入CT或MRI扫描的数据, 进行图像的分割, 阈值的设定图像增强等处理, 直接建立三维模型, 再使用Geomagic逆向工程软件, 对模型进行修饰、剪裁、光顺等后处理操作, 使模型达到科研要求, 再导入到有限元分析软件, 如ANSYS、ABAQUS, 进行相关的数据分析。此方法比CT建模简单, 省略了用图形处理软件对数据进行边界提取、二次处理等大量的工作[19]。

1.6 数字化虚拟人的建模方法 数字化虚拟人是现代计算机信息技术与医学等学科相互结合的成果[20]。以人体组织连续断面构成数据集合, 可以逼真地建立骨骼、肌肉、血管等器官组织的数字化模型。

2关节盘及关节窝的模拟

在模拟介于关节窝和髁突之间的关节盘时, 国内外学者所模拟的厚度有所不同。大部分学者模拟的关节盘的厚度为2mm, 数据源自Hansson[21]所做的研究。Nagahara等[14]所建的模型中, 关节盘前、中、后带均为2mm;Tanaka等[4]建模研究紧咬牙时, 颞下颌关节不同区域的应力分布时, 所模拟关节盘厚度为2mm;Cheng等[22]研究颞下颌关节紊乱的患者时所见模型关节盘厚度为2mm;Groning等[23]建立的模型, 颞下颌关节软组织的厚度选择为3mm;胡凯等[7]建模时, 关节盘厚度也为2mm;安虹等[24]在构建功能状态下无牙颌下颌骨及颞下颌关节三维有限元模型时, 关节盘厚度也为2mm;周祺[25]在研究下颌升支矢状劈开后退术对颞下颌关节影响时建立的颞下颌关节也选择平均厚约2mm的关节盘以及长方体结构模拟颞下颌关节窝的皮质骨、松质骨等部分, 从而建立完整的颞下颌关节的三维有限元模型;刘展等[13]认为关节盘位于髁突和颞骨之间, 除了具有关节软骨层的作用外, 还能保持与髁突、颞骨的接触, 协调关节的运动, 是维持TMJ稳定的重要因素, 根据其解剖特征和关节表面形状建立关节盘的实体模型, 前、中、后带的厚度分别为2、1、2.7mm;周学军等[26]建模时, 模拟的关节盘中带1.2mm, 关节盘后带2.0mm, 于关节盘表面建立带凹的方块皮质结构代表TMJ的关节窝结构;胡林华等[10]建立的颞下颌关节-下颌骨-Herbst矫治器系统三维有限元模型, 模拟的关节盘厚度中带、后带分别为:1.2、2.0mm。

3肌肉及韧带的模拟

史真等[28]建立下颌骨牵张成骨三维有限元模型, 对咀嚼肌、下颌韧带采用杆单元模拟其约束, 杆单元材料定义为只受拉不受压的非线性材料, 单元横截面积与各自模拟的肌肉和韧带截面积相同;周学军等[11]模拟咀嚼肌被动受张作用下颌骨时, 采用缆索元即无间隙的受拉单元于附着中心模拟肌肉、韧带的约束, 更接近实体状况。

4关节盘与髁突、颞骨之间的模拟

三维有限元法 篇5

根据航空等领域内锪窝铆接及锪窝锣接构件的典型结构特征,采用三维的十节点四面体等参有限单元模型,分别对无裂纹及孔边含裂纹锪窝孔/直通孔结构进行了模拟分析;得到了锪窝孔构件的危险部位及90°, 120°锪窝孔边扇形角裂纹的应力强度因子,给出了覆盖面广的计算曲线; 通过对计算结果的分析,讨论了裂纹长度、孔径以及板厚等因素对应力强度因子的.影响.和已有的文献比较表明,本文数值结果精确,方法可靠.

作 者：沈海军 郭万林  作者单位：南京航空航天大学,飞机系,江苏,南京,210016 刊 名：航空学报  ISTIC EI PKU英文刊名：ACTA AERONAUTICA ET ASTRONAUTICA SINICA 年，卷(期)：2002 23(2) 分类号：V215.6 关键词：锪窝孔   有限元   应力强度因子   扇形角裂纹  

★ 建筑材料的强度和比强度是什么？

★ 电场强度教学设计

★ 孔和尚 语录

★ 孔老师作文

★ 孔林导游词

★ 低压电电工工转正申请书范文

★ 海州湾特征因子与赤潮形成初探

★ 孔学堂参观记作文

★ 孔祥熙简介

四两拨千斤 铅球三维训练法 篇6

多数教练员、考生专注于对技术要领的教授和学习，忽略了对铅球运行轨迹和用力三维方向的全面认识。而高中生已经具备较好的三维空间认知能力，如在理论和训练课中引入推铅球用力的三维方向，训练就能有的放矢，学生将有限的力量使用在正确的方向和轨迹上，发力就能起到四两拨千斤的高效，从而大幅提高原地推铅球的成绩。

[关键词]高考；推铅球；用力；三维；方向

一、原地推铅球用力的三维方向

如图1（左、右脚着地点和铅球投影点的平面图）：以右手推铅球为例，图中箭头为投掷方向，1、2分别为左、右脚，3为推球者所持铅球投影点。图中左右脚与铅球投影点连接，形成稳固的三角形。

如图2（铅球场地和三维坐标侧视图）：以运动员直立，面向落地区出手方向时的髋关节中心点建立三维坐标原点O，铅球水平出手方向为X轴，垂直于地平面为Y轴，垂直于XY平面为Z轴。铅球运行轨迹分别投影于XY平面（如图3）和XZ平面（如图4）。人体作用铅球的距离AB用实线表示，铅球出手后的飞行轨迹BC用虚线表示。

我们会发现铅球在人体作用下作螺旋上升、离心运动。将铅球运行轨迹分别分解为XYZ三个方向，运行轨迹是：后→前、下→上、右→左的三个方向。为便于直观地表述和教学应用，在此我们将其称为铅球运行轨迹的前后方向、上下方向和左右方向。

（一）原地推铅球的前后方向

最大远度是每名考生追求的目标，学生往往盲目关注远度，毫无方向地出手，多数抛物线过低。学生仅仅关注向前，而忽略由后向前的过程，也就是前后方向。

如图3所示：人体作用铅球的距离为AB， AB投影在X轴上为AX1，投影在Y轴上为OY1。AX1即为原地推铅球的前后方向。推铅球出手最佳角度是38-42度，则∠OAB的角度也为38-42度，根据三角形的边角关系可以判定，AX1>OY1，说明人体作用在铅球上水平距离上大于垂直方向上的距离。要使铅球达到最大远度，首先要使铅球在前后（水平）方向上达到最大的位移距离，其次才是上下（垂直）方向。

怎么才能获得最大前后位移呢？简化动作，弱化屈膝蹬地，降低发力难度，强化前后方向，充分感受前后位移。1.站位适度拉开：左脚与右脚的前后距离适度拉开。以身高170CM学生为例，左右脚前后距离为80--90CM为最佳。 2.铅球起始位置远离落地区：在“原地”规则的限制下，能最大限度增加作用距离的方式就尽可能地使铅球投影点接近后半圆，远离投掷方向。使用原地背向、深蹲右脚的姿势，能获得最大的前后位移。3.铅球出手点接近落地区：铅球出手瞬间的投影点应该在抵趾板正上方，或者是落地区的扇形内，使铅球在最大限度的前后方向的轨迹内运行。

（二）原地推铅球的上下方向

如图3所示：人体作用铅球的距离为AB， AB投影在Y轴上的距离为OY1。而OY1为人体作用在铅球的上下方向，OY1即为人体作用在铅球上的垂直位移。Y1与地面的距离即为铅球出手的高度，即身高+臂长=铅球出手高度。

怎么才能获得最大的上下位移呢？简化动作，降低难度，弱化前后，强化上下，强化蹬地，强化出手速度和高度，充分感受上下位移。1.橡皮筋阻力练习法：将皮筋的一端固定在学生后下方，学生右手持皮筋的另一端进行徒手模仿练习。2.平行站位推高练习法：持球正对投掷方向，两脚平行站位，下蹲起，体会蹬地发力由下向上的运行轨迹。3.侧向屈体“半屈膝”推高练习法：持球侧对投掷方向，左腿略屈，右脚半屈膝，弯腰含胸，将铅球推至一定的高度。

（三）原地推铅球的左右方向

如图4所示：为人体作用铅球的距离AB投影在XZ水平面上的轨迹，不难发现其起始阶段是一条弧线。弧线是由于推铅球技术中的“转”体动作和上体向前上方的运动综合所产生。往往这又是在日常训练中最容易忽视的一个方向，却又是最大限度增加作用距离和提高铅球出手初速度最重要的因素。

怎样才能获得最大的左右位移呢？简化动作，弱化前后和上下方向，强化左右方向，强调肩轴的旋转过渡，强调初速度的左右方向，充分感受右肩的左右位移（旋转）。1.平行站位半蹲转体练习法：持球正对投掷方向，兩脚平行站位，半蹲同时右转体约45度，强调肩轴的旋转，利用蹬地和向左转体45度的合力，将铅球掷出。2.前后站位半蹲转体练习法：持球正对投掷方向，两脚前后站位，右脚尖与左脚后跟齐平，半蹲同时右转体约45度，强调肩轴的旋转，利用蹬地和向左转体的合力，将铅球掷出。3.平推铅球：弱化上下方向和出手角度，利用前后方向的初速度，强化右肩由右向左的旋转。

二、三维用力方向的实际运用

简单地说，铅球完整运行方向就是以上三个方向的合成，这三个方向缺一不可，相辅相成。开始姿势：右手持球于锁骨窝，左手侧上举，抬头挺胸顶髋；左脚直立，左脚脚尖和身体正对投掷方向，右脚前脚掌着地；身体重心落在左脚。身体正对投掷方向的预摆，利于学生不断重复地体会发力顺序和技术，在抬头挺胸的高位发力出手，形成正确的出手高度。重视预摆和左侧摆停的技术。“预摆”就是开始姿势→准备姿势→开始姿势的循环过程。在预摆中，包含有推铅球下肢发力核心动作“蹬、转、挺”，预摆使人体和铅球处于最低位和高位，是人体用力作用铅球的初始轨迹，运行路线长；左手侧上举和身体对正摆停的左侧技术，利于提高出手速度和增加做功距离。

参考文献：

[1]引自期刊：作者：李红莲、陈建忠《原地背向“旋展式”推铅球技术介绍》，《体育教学》，2008年第6期，第53页.

[2]引自期刊：作者：瞿昂《我国优秀男子大学生铅球运动员左侧技术的三维运动学分析》，《浙江体育科学》，2008年第5期，第123页.

三维有限元法 篇7

1 资料与方法

1.1 一般资料:

以我院在2010年3月至2013年7月收治的68例口腔正畸矫治患者为研究对象, 随机将患者分为两组, 每组各34例病例, 观察组采用三维有限元法治疗, 男性19例, 女性15例, 年龄在17~56岁, 平均年龄 (32.26±9.75) 岁。对照组采用常规治疗法治疗, 男性18例, 女性16例, 年龄在18~55岁, 平均年龄 (33.75±9.14) 岁。两组患者在一般资料上无统计学意义 (P>0.05) 。

1.2 治疗方法

1.2.1 对照组:

以患者畸形状况为依据进行矫治, 患者需行X线片检查, 并做好牙模准备工作, 同时还需实施分牙、拔牙等各项操作, 完成上述操作后, 便可将戴上牙套, 便于矫正牙齿, 待牙齿处于干燥状态后, 可蚀刻, 并将粘结材料加入其中, 选取弓丝, 结扎于牙套上, 完成矫正30 d后, 需返院复诊, 对其定型与调整后, 可摘除牙套。

1.2.2 观察组:

利用有限三元法治疗, 首先需做好网格划分操作, 并对待解区进行分割, 确保其属于有限单元集合, 网格的主要由四面体、节点、单元组成, 分析节点力与节点位间的关系, 分析多单元形成整体, 并明确节点位移数据。建数学模型, 并使用CT仪对二维图像进行获取, 并将图像数据, 利用三维有限元软件进行分析后, 便可形成三维有限模型。以患者上下颌骨、牙体组织、颅颌面修复情况为依据, 对矫治方式进行选择, 主要利用mbt直丝弓矫治器, 完成矫治时间达半个月后, 需返院复诊。

1.3 矫正效果评价标准。

优:患者完成矫正治疗后, 对口腔部位与面部均起到了美化作用, 无畸形现象, 且未出现不良反应, 美化效果非常好。良好:患者经矫正治疗后, 无不适度感, 畸形现象有明显改善。差:患者接受正畸治疗后, 矫正未起到明显作用, 与治疗前相较变化不大。

1.4 统计学方法:

对本组研究的数据采用SPSS16.0统计软件进行分析, 计量资料以均数±标准差 (±s) 表示, 采用t检验, 对计数资料采用χ2检验, P<0.05为差异有统计学意义。

2 结果

见表1。观察组中有26例矫正情况为优, 7例患者矫正情况为良, 1例患者矫正情况较差, 总优良率为97.1%;对照组中有20例矫正情况为优, 8例患者矫正情况为良, 6例患者矫正情况较差, 总优良率为82.4%。观察组矫正优良率高于对照组, 对比差异显著, 具有统计学意义 (P<0.05) 。

3 讨论

有限元分析法目前在口腔生物力学中的应用非常广泛, 这种分析法于1954年被提出, 关键在于建立模型, 模型的力学、几何相似性与网格划分, 均关系到计算结果。口腔组织的组成非常复杂, 且形态具有不规则的特征, 在建模与采集数据方面难度较大。国内与国外学者都致力于研究出一种有效的建模方法, 近年来, 三维有限元法在口腔建模中有着良好的应用, 并且拓展了应用范围, 利用它可建立逼真的牙槽骨牙周模型与牙体组织模型[4]。

有限元力学模型的基础就是实验, 首先需要对牙颌组织进行获取, 并以牙颌组织为依据, 建立组织模型, 在获取信息的过程中, 需尽量准确, 确保三维空间几何数据的准确性, 这对于建模具有重要意义。国外有多种获取方式, 其中主要包括利用切片、磨片获取, 或者通过非接触式三维形态获取。牙体组织的计算与模拟较为复杂, 将三维有限元法应用于口腔正畸矫治中, 可有效改善患者病情, 这对于患者疾病的好转具有重要意义[5]。

畸形矫治的主要原理为将关外力施加于畸形颌骨上, 或者将异常肌力去除, 在牙周组织、机体颌骨等组织内部会有生物力学反应产生, 便于改建产生组织学, 可确保牙颌系统的正常发育以及功能平衡性。随着口腔技术的进步与发展, 种植体被应用于口腔正畸矫治中, 使用三维有限元法对患者的腭部支抗种植体进行研究, 当处于受力状态时, 支抗种植体会做复合运动。据相关资料显示, 即刻垂直加载的过程中, 种植体周围无法将骨应力集中, 且无法均匀分布, 于颊舌部位加载的同时, 其应变、应力会有增大的迹象, 且呈现为分布不均的情况[6]。现阶段, 关于周围骨组织与种植体如何结合, 还有待进一步研究。

有限元法是一种将理论与现代信息技术相结合的分析方式, 在口腔生物力学中得到了有效应用, 可前景十分开阔。因人体组织具有材料非线性、结构不规则等特征, 对计算机制图能力与其性能也提出了较高要求。现阶段, 在口腔领域中, 有限元建模技术还有很大的提升空间, 它已经由静态转变为动态, 具备很多优越性, 以后还将有更好的发展。若能够将三维有限元法同其他方法相结合, 便能够使研究更为深入。

从本次研究中可看出, 观察组患者采用三维有限元法治疗, 对照组采用常规法治疗, 观察组中有26例矫正情况为优, 7例患者矫正情况为良, 总优良率为97.1%, 这表明采用三维有限元法对患者进行矫治, 具有较高的安全性与可靠性, 便于了解患者口腔内的具体情况, 通过建模将患者口腔状况变得更为明确, 有针对性的采取治疗措施, 可大大提高矫治优良率, 患者对于这种方式也更乐于接受, 值得临床推广应用。

摘要：目的 研究三维有限元法在口腔正畸矫治中的应用价值。方法 选取我院在2010年3月至2013年7月收治的68例口腔正畸矫治患者的临床资料, 随机将患者分为两组, 每组各34例, 观察组采用三维有限元法治疗, 对照组采用常规治疗, 比较两组患者的治疗效果。结果 观察组中矫正情况优者26例, 矫正情况良者7例, 矫正情况差者1例, 总优良率为97.1%;对照组中矫正情况优者20例, 矫正情况良者8例, 矫正情况差者6例, 总优良率为82.4%。观察组治疗优良率优于对照组, 对比差异显著, 具有统计学意义 (P<0.05) 。结论采用三维有限元法治疗, 可取得更为满意的矫正效果, 患者更易于接受这种治疗方式, 值得临床推广应用。

关键词：三维有限元法,口腔矫正,临床疗效

参考文献

[1]赵玺, 韩蕊, 米丛波.三维有限元法在口腔正畸生物力学研究中的应用[J].中国实用口腔科杂志, 2011, 4 (2) :116-118.

[2]李冬梅, 李金源, 冯宝华, 等.牙颌组织三维有限元模型在口腔正畸中的研究进展[J].当代医学, 2011, 16 (7) :24-25.

[3]武倩倩, 赵桂芝, 柯杰.三维有限元法在口腔正畸学中的应用进展[J].空军医学杂志, 2012, 28 (1) :24-25+45.

[4]吴志芳, 雷勇华.三维有限元法在口腔正畸中的应用和展望[J].国际口腔医学杂志, 2013, 40 (6) :804-807.

[5]张莹, 彭惠.三维有限元法在正畸生物力学中的研究进展[J].齐齐哈尔医学院学报, 2013, 34 (2) :249-251.

三维有限元法 篇8

由于管道消声系统内部存在气流,且气流的存在对管道和消声器内部的声传播有着直接影响,所以在消声器声学性能计算及分析中有必要考虑气流的影响[1,2,3,4]。由于实际消声器的结构比较复杂,其内部流场和声场本质上是三维的,为了精确预测其声学特性需要使用三维数值方法。有限元法是一种先进高效的数值计算方法,已被广泛应用于求解各种类型的工程问题。然而声学有限元法目前还只能用于求解无流[5,6]、具有低马赫数非均匀流[7]和均匀流[8]时管道及消声器的声学问题。对于具有较高马赫数亚音速流时管道及消声器声学问题计算的有限元法目前尚无研究成果发表。本文的主要目标就是将有限元法应用于计算三维势流中的声传播问题,推导声学有限元法计算公式并给出数值求解过程。作为应用实例,将三维有限元法应用于计算和分析具有三维势流时管道和消声器的声学特性。

1 流场与声场的控制方程

内燃机排气消声器内的气流马赫数基本在0.3以下,可视为不可压缩有势流动,于是介质的流动速度可表示为:

式中,φF为流速度势。

流体连续性微分方程为:

式中,ρ为介质密度。对于不可压缩流体,ρ不随时间和空间位置而变化,于是有:

将式(1)代入式(3)可得到流场的控制方程,即Laplace方程:

在三维势流中,稳态声场的控制方程[9]为:

式中,A为声速度势;M为气流马赫数,M=V/c0,c0为声速;γ为比热比;j为虚数单位,。将M=V/c0代入式(5)的第五项和第六项,考虑到式(3)成立,并且忽略掉马赫数的三阶以上高阶小量,则式(5)可简化为:

将式(6)展开后得:

为了表示方便,引入以下记法:

则式(7)可简记为:

式中,Mx、My、Mz分别为马赫数M在x、y、z方向的分量;I13={1 1 1}是一个1×3的单位阵。

声压p和质点振速u与声速度势φA之间的关系为:

式中,ρ0为介质密度;ω为角频率。

由于稳态声场的控制方程中含有马赫数,所以为获得三维势流中的声场分布,需要首先求解方程(4)得到流场信息,然后将结果代入方程(8)求出声场信息。

2 流场计算的有限元法

在流场的有限元计算中,选择一系列的插值函数来表示稳态的流速度势:

式中,{φF}e为节点流速度势列向量;{N}为形函数列向量;上角标T代表列向量的转置。

将式(11)代入式(4),使用迦辽金加权残数法可以得到流场的有限元方程为:

消声器内流场计算中的边界条件为:

(a)消声器壁面上法向流速为零

(b)消声器进口设为速度边界条件

(c)消声器出口设为速度势边界条件

将边界条件式(13)、(14)代入方程(12)可得:

其中,

将系统中的所有单元进行组装后可得:

式中,[M]、[F]分别由相应的单元系数矩阵扩充而成;{φF}由所有节点处的流速度势构成的列向量。应用边界条件(15),求解方程(19)可得到各节点处的流速度势。为了求出各节点处的流速,将式(11)代入到式(1),于是有:

3 声场计算的有限元法

在声场有限元计算中,选择插值函数来表示声速度势,即:

式中,{φA}e为节点流速度势向量;{N}为形函数列向量。

将式(21)代入式(8),使用迦辽金加权残数法可得:

引入以下记法:

式中,边界表面Se=Ser+Sev+Sez,Ser、Sev和Sez分别代表刚性壁面、已知法向质点振速的边界表面和已知法向声阻抗的边界表面。则式(22)可简记为:

上式即为声场的有限元方程。

消声器声场计算中的边界条件为:

(a)消声器壁面为刚性,法向速度为零

由此可得:

(b)消声器进口设为质点振速边界条件,本文设为un=1,则有:

由此可得:

(c)消声器出口设为无反射端,即:

将式(32)代入式(26)可得:

其中,

将(a)、(b)和(c)3个边界条件代入式(26)可得:

将式(35)代入式(27),整理可得:

将系统中的所有单元进行组装后得到:

式中,[KA]、[MA]、[DA]、[Cz]和{Fin}分别由相应的各单元系数矩阵扩充而成;{φA}为所有节点处的声速度势构成的列向量。值得注意的是:向量{Fin}仅由消声器进口面贡献,[Cz]仅由消声器出口面贡献,求解方程(37)即可得到各节点处的声速度势。联合已求得的流场变量,通过式(10)即可计算出各节点处的声压值。

4 传递损失的计算

传递损失定义为消声器进口处的入射声功率级与出口处的透射声功率级的差值。由于传递损失只与消声器本身结构有关,不受声源特性和尾管辐射特性的影响,因而是消声器研究中最常用的声学性能指标。消声器进出口管的直径通常较小,截至频率较高。当进出口管道内满足平面波条件时,消声器的传递损失[1]可表示为:

式中,S1、S2分别为消声器进、出口的横截面积;pi和pt分别为消声器进口处的入射声压和出口处的透射声压。由于出口为无反射端,只包含前进波pt,而在进口处同时存在入射波pi和反射波pr,所以消声器的传递损失可表示为:

式中,p1、v1分别为消声器进口处的声压和质点振速;p2为消声器出口处的声压。当消声器进口处的质点振速v1给定时,使用有限元法可计算出进口和出口处的声压p1和p2,于是代人式(39)即可求出消声器的传递损失。

5 计算实例与分析

为了计算有流时管道和消声器的声学特性,首先需要计算流场分布,然后计算声场分布。根据进出口处声压和质点振速的计算结果即可求出相应的声学特性。在本文中,将计算区域划分成10节点四面体二次单元,在进出口处由三维单元自由面生成相应的二维单元,本文为6节点三角形二次单元。取空气的密度ρ0=1.2kg/m3,声速c0=340m/s。

5.1 直管道的四极参数

为了验证本文方法的正确性,首先需要计算圆形截面直管道的四极参数。本文取管道长度为0.4m,直径为0.1 m,计算得出管道的截止频率为4148 Hz。在平面波范围内,可以通过有限元计算结果与一维解析解的比较来验证本文方法的正确性。图1为直管道的四极参数。图1a比较了无流时使用传统声学有限元法和一维解析方法计算得到的直管道四极参数。图1b比较了进口处气流马赫数为0.3时使用本文有限元法和一维解析方法[1]计算得到的直管道四极参数。结果表明:二者完全吻合,本文方法预测具有较高马赫数流时管道声学特性是正确的。比较图1a和1b可以看出介质流动对声学特性的影响,由于本文方法考虑了介质流动效应,因而更加真实地描述了流动介质中声传播的物理现象。

5.2 消声器的传递损失

使用本文有限元法计算如图2所示(单位为mm)的具有外连管的双级膨胀腔消声器的传递损失,考察三维势流对消声器声学性能的影响。不同气流马赫数时消声器传递损失的有限元计算结果如图3所示。介质流动对消声器传递损失的影响主要体现在运流效应。由于介质的流动使得消声器内的入射声波与反射声波的传播速度不同,因而影响了消声器的消声频率特性。比较图3中M=0和M=0.3时的传递损失,可以看出明显的频移现象,频率越高、流速越大,这种频移越明显。同时还可看出,随着气流马赫数的增加,双级膨胀腔消声器的传递损失曲线明显地向低频方向移动,因而改变了各个频率下的消声量。

为了进一步考察气流对复杂结构消声器声学性能的影响,使用本文有限元法计算并分析如图4所示(单位为mm)的具有两个内插管的双级膨胀腔消声器的传递损失。消声器的传递损失随气流马赫数的变化规律如图5所示。介质流动对该消声器低频声学特性的影响较小,对高频特性的影响较大。传递损失曲线随进口气流马赫数的增大而明显地向低频移动,进而改变了各频率下的消声量。因此三维势流对复杂结构消声器声学性能的影响是应该加以考虑的。

6 结论

(1)建立了有三维势流存在时管道和消声器声学特性预测的有限元法,并通过实例验证了该方法的正确性和计算精度。

(2)与传统的声学有限元法相比,新方法考虑了控制方程中气流马赫数的二阶小量,因此适用于求解具有较高马赫数亚音速流中的声传播问题。

(3)对双级膨胀腔消声器的传递损失计算结果表明:随着气流马赫数的增加,消声器传递损失曲线向低频方向移动,因而改变了各频率下的消声量。

(4)由于有限元法不受几何形状的限制,且能够考虑消声器内的三维流场和三维声场对消声效果的影响,因而能够用于计算不规则形状或复杂结构消声器的声学性能。

摘要：将有限元法应用于预测有三维势流时管道和消声器的声学特性,推导了三维势流中声传播问题的有限元法计算公式,给出了管道和消声器声学特性的计算方法和数值实施过程。与传统的声学有限元法相比,新方法考虑了声学控制方程中气流马赫数二阶小量的影响,因此适用于求解具有较高马赫数亚音速流中的声传播问题。使用三维有限元法首先计算了管道的四极参数,并与一维解析解比较验证了新方法的正确性和精度,然后使用三维有限元法计算双级膨胀腔消声器的传递损失。数值计算结果表明:随着气流马赫数的增加,消声器传递损失曲线向低频方向移动,因而改变了各频率下的消声量。为了精确预测消声器的声学特性,三维势流对管道和消声器内声传播的影响应该加以考虑。

关键词：内燃机,消声器,声学特性,有限元法,三维势流

参考文献

[1]Munjal M L.Acoustics of ducts and mufflers[M].New York:Wiley-Inter-Science,1987.

[2]Ji Z L,Ma Q,Zhang Z H.Application of the boundary ele-ment method to predicting acoustic performance of expansionchamber mufflers with mean low[J].Journal of Sound and Vi-bration,1994,173(1):57-71.

[3]Ji Z L,Ma Q,Zhang Z H.A boundary element scheme for e-valuation of four-pole parameters of ducts and mufflers with lowMach number non-uniform flow[J].Journal of Sound and Vi-bration,1995,185(1):107-117.

[4]王雪仁,季振林.有流消声器四极参数和传递损失的边界元法预测[J].计算物理,2007,24(6):717-724.Wang X R,Ji Z L.Prediction of four-pole parameters andtransmission loss of silencers with flow by the boundary ele-ment method[J].Chinese Journal of Computational Physics,2007,24(6):717-724.

[5]Young C I J,Crocker M J.Prediction of transmission loss inmufflers by the finite element method[J].Journal of the A-coustical Society of America,1975,57(1):144-148.

[6]Young C I J,Crocker M J.Acoustic analysis,testing and de-sign of flow-reversing muffler chambers[J].Journal of the A-coustical Society of America,1976,60(5):1111-1118.

[7]Peat K S.Evaluation of four-pole parameters for ducts withflow by the finite element method[J].Journal of Sound and Vi-bration,1982,84(3):389-395.

[8]Tsuji T,Tsuchiya T,Kagawa Y.Finite element and boundaryelement modeling for the acoustic wave transmission in meanflow medium[J].Journal of Sound and Vibration,2002,255(5):849-866.

深大基坑开挖三维有限元分析 篇9

深基坑支护设计最根本的要求是“安全、造价低、施工方便”,这是全面考察主体结构的表现、水位高低的波动、现场环境的变动等,对多种支护结构设计进行对比,从中选择最适合该施工情况的支护形式。深基坑工程是危险性工程,稍有不慎,会对人们的生命安全和财产安全造成不可挽回的影响,甚至造成不可估量的后果。因此,在开发利用地下空间时,需充分考虑基坑开挖的安全,尤其是深大基坑的施工。在深基坑开挖过程中如何控制基坑开挖的位移和稳定性是其施工难点。

1 工程概况

某深大基坑工程位于市区,该项目拟建3栋商住楼,其中塔楼23层,裙楼5层,地下3层。该基坑的土层以粉土和砂土居多,并且水位常年很高且变化不明显,地质土体不好,土的弹性模量很高,强度低,容易产生变形;砂土透水性好,但稳定性差,边坡坍塌、涌砂、涌水现象在施工过程中容易产生,围护结构没有大的刚度不能实现支护效果,止水效果不明显,难以进行施工。考虑水文地质条件和经济成本原因,结合专家意见,选取钻孔灌注桩方案,桩间距为1 250 mm,外加三轴搅拌桩,采用850 mm止水帷幕。为了防止基坑开挖引起的钻孔灌注桩侧向位移过大、地基土变形,合理确定基坑的开挖工况、支撑的程序尤为重要。深基坑的开挖深度与最下道钢支撑到坑底的距离不能过大。在全面整合研讨之后,可以采取三步来竖向开挖深基坑,配置3道钢支撑,支撑材料选取直径106 mm、壁厚16 mm的钢管。

2 深基坑支护设计

2.1 深基坑支护原理分析

随着计算机的发展与普及,有限单元法以其独有的优势应用于挡土结构分析,它易于模拟、灵活、约束少,还能有效地考虑基坑开挖过程中的多种影响因素。

有限单元法也有其明显的不足,如果模型选取不当,或者简化的模型与现实差别太大,以及模拟参数的选取没有足够的准确性,其计算结果的可信度就值得怀疑。因此,要想得到准确的模拟结果,就必须反复地对比模型的选取是否合适,模型简化是否与现实接近,材料的本构模型以及参数选取是否真实反映实际情况等。

2.2 深基坑支护计算模式

利用有限元方法来计算板桩结构的原理及模式如图1所示,其主要原理是将弹性地基梁绕中心位置旋转90°,然后将挖面以上的土体弹簧单元移动到非工作结构处。墙后具有支撑墙体的结构,该结构在支撑的过程中同时受到主动土压力的作用,对于墙体前需要施工的土及钢结构的支撑结构为弹性地基。根据就近单元原则对墙体上的土承受的压力进行排列划分,从而可将单元简化为单元节点上的一点,然后对该点的集中力载荷进行分析。对于同样区域的土下部分土体的弹簧载荷分析的原理与墙体类似,也是对土体的就近单元进行划分,然后将载荷分析简化为对某个节点的集中力荷载分析,其节点上承受弹性支撑力,然后利用有限元的相关原理和计算方法对支护墙体进行载荷研究和稳定性分析。上述的分析过程可以简要地表述如下:首先将分析土域的土体以及其支撑架构当作弹性支撑,然后利用墙体与支撑结构的共同协调作用可知它们将同时变形,因此利用同时变形的原理对支护结构进行分析的结果是相对比较准确的。

采用排桩加内支撑的结构形式确定基坑开挖施工工况是确保工程进行的前提,应根据实际施工过程中基坑开挖的进度情况,来确定支护结构在开挖过程中的受力计算简图。

在进行支护结构内力分析时,应针对各个不同的施工阶段的工况进行全面的分析计算,才能确保支护结构在全施工过程中的安全性。

3 深基坑支护有限分析

土体是由固体、液体、气体三相结构组成的一种复杂体,土体的本构关系主要是指其对应的应力与应变之间的关系。土体的组成元素相对比较复杂同时成分之间的离散性较大,由于元素之间的物理性能差异较大,因此对于土体整体来讲具有非线性、粘塑性等物理性质。国内外对土体的本构关系研究已经有几十年的经验,在相关的研究理论的发展中,专家提出了上百种关于土体本构的研究模型,同时在实际的工程应用实践中已经将成熟的研究理论引入其中,同时利用先进的计算机编程的方式来建立相关的数学模型,逐步引入数值分析技术以解决一些工程中的实际问题。在土体相关关系从理论研究发展到实际应用的过程中,成熟的理论及研究模型要被人们认可并且应用于实践中是十分有限的。

3.1 深基坑开挖有限元模拟过程

该施工模拟过程主要分三步进行,其在开挖过程中主要采用“空单元”法进行。其主要的计算过程及步骤如下:

对深基坑在施加足够重力载荷的同时,计算其初始的应力场和位移场。

第一次开挖深度为1.2 m,开挖过后的深基坑单元用“空单元”来表示,并对该过程中产生的应力场和位移场进行分析和计算。

设置第一道支撑,第二次开挖深度为6.7 m,开挖过后的深基坑单元用“空单元”来表示,并对该过程中产生的应力场和位移场进行分析和计算。

设置第二道支撑,第三次开挖深度为12.2 m,开挖过后的深基坑单元用“空单元”来表示,并对该过程中产生的应力场和位移场进行分析和计算。

设置第三道支撑,第四次开挖深度为15.5 m,开挖过后的深基坑单元用“空单元”来表示,并对该过程中产生的应力场和位移场进行分析和计算。

深基坑支护有限元分析模型如图2所示。

3.2 深基坑开挖有限元分析结果

深基坑支护结构水平位移对比图见图3。

对图3的位移图进行分析可知,围护桩支护体系的水平位移在一定的范围内,本次模拟试验中其水平位移在30 m以内,其中在第一次开挖过程中基本不发生水平位移,第二次开挖过程中可能在坑底处产生8 mm左右的水平位移,第三次开挖出现的水平位移高达18 mm,本次水平位移发生的位置集中在距离坑深6 m左右的位置,而第四次开挖过程中产生的水平位移达到23 mm左右,该水平位移主要发生在距离坑底8 m处。

在各次开挖过程中,围护桩的弯矩会发生如图4所示的变化。

在具体的工程实践中,例如本次的深基坑开挖工程中围护桩承受土体的水平载荷而发生的弯矩会随着开挖深度的加深不断增大,同时支撑土体的支撑结构受到的弯矩也在不断增加。通过以上分析可知,最大弯矩值发生在第四次开挖过程中,其值为-1 839 k N·m,将该弯矩值与围护桩所能承受的弯矩极限值比较,整个开挖过程中围护桩的结构性能基本上能得到保障。同时在施工过程中随着开挖深度的不断增加,其支撑的支护桩的正、负向的弯矩也会不断地增大,因此为了保障施工过程的安全性,必须增加施工的截面尺寸,同时向支护桩中配置钢筋以增强其支撑力,另外围护桩中的双向配筋也能起到对施工过程中的安全性能增强的作用,能够充分抵抗桩身弯矩的变化而产生的影响。

对于支护桩的安装也需要满足一定的安全要求,例如支护结构应该埋置的深度的问题,如果埋深增加,支护桩的侧向位移将会减少。当支护桩的埋深得到保障后,在施工过程中其发生的最大位移一般处于开挖面附近,如果此时继续对支护结构增加埋置深度,其抵抗基坑变形的能力并不会得到明显的提高。国内外相关专家的研究发现,在保证基坑底部不隆起以及支护桩的埋置深度满足条件最小的情况下,其围护桩的位移量与支护桩的埋置深度并没有太大的关系。因此,在选择支护桩的埋置深度时主要考虑其正、负向的弯矩值,以保证在施工过程中能够充分利用材料实现对桩的优化设计。

在开挖过程中深基坑底部出现隆起的主要原因是在土体卸载过程中其原始的应力状态发生了相应的变化。在整个开挖过程中其隆起值的曲线图如图5所示。由图5可知,其基坑中心为水平距离的零点,基坑的总宽度为18 m。

对图5的曲线图进行分析可知,基坑底部的隆起量会随着开挖深度的增加而不断增加,大约在基坑深度的中间位置其发生的隆起量最大,基坑侧壁依靠围护桩的保护以及约束,因此其隆起量相对较小可以忽略。第四次开挖过程中基坑深度的中间位置的隆起量为7.86mm,本次的工程实践中坑底隆起量的变化基本上是均匀的,没有出现较大幅度的隆起起伏。因此,整个施工过程中坑底隆起处于弹性工作状态,整个工程环境安全,不会因为隆起而发生基坑破坏。施工中围护桩外部的土体的变形曲线如图6所示。

4 结语

在水平压力的作用下,开挖后的基坑会使支撑的围护桩发生一定的变形,从而使得整个施工区域除了围护桩以内的土体出现沉降。随着开挖深度的加深,土体的沉降量越来越大,其最大的沉降量为21 mm,而桩附近的土体因为受到围护桩或者支护桩的保护其沉降量不明显;施工土体距离基坑越远,开挖对其沉降的影响也越小,其在坑深3倍处位置的变形量趋于平稳。

水是影响基坑安全性最重要的因素。本文在建立有限元分析模型时对水位进行了设置,研究了水位变化的情况下基坑围护结构的位移和基坑外侧土体的沉降,可供类似工程参考。

摘要：深基坑支护在城市建设中使用广泛,如何保证深基坑的安全性、合理性已经成为深基坑设计中的重要问题。结合工程实例,利用岩土计算软件理正深基坑7.0和有限元分析软件Midas-GTS对基坑的安全性进行了分析,其分析结果可为优化设计和施工提供参考,为类似工程提供借鉴。

三维随机渗流摄动有限元模拟研究 篇10

随着近几十年来计算技术及可靠性理论的发展,岩土工程数值分析中不确定性的考虑越来越受到重视。岩土渗流计算作为岩土工程数值计算分析的一个重要分支,存在很多不确定因素,因此,将渗透参数的随机性引入到渗流计算中是十分必要的。现有的对渗流场进行随机分析中Monte Carlo法[1]是比较常用的方法,该方法需对渗透参数按照一定分布形式进行大量的随机抽样,通过对确定性渗流分析的多次重复求解,从而得到渗流场的统计特征。此方法计算量比较大,限制了其在实际工程中的应用。

摄动法是用于分析非线性问题的强有力的方法,由于理论明确、计算量小,摄动有限元法在力学及其他工程科学领域有广泛的应用。80年代,小参数摄动技术被引用到有限元计算中[2,3],Erik Vanmarcke[4]提出随机场的局部平均理论并将它引入随机有限元法,使得摄动随机有限元法更趋于实用。90年代,陈虬与刘先斌[5]编制了计算平面弹性问题的二维摄动有限元法程序(PS FE2 D E)。渗流问题与弹性问题有很大的相似性,本文参考原程序,将渗透张量视为各向异性的平稳随机场,各介质呈非均匀分层分布,应用随机场的局部平均理论对渗流场进行离散,得到一系列随机变量来近似模拟空间渗流场的随机性。最后,编制相应的三维渗流计算的摄动随机有限元程序,计算随机渗流场,将此计算结果与蒙特卡罗法计算的渗流场进行比较,验证此方法的正确性和可行性。

2 三维渗流随机场的局部平均离散

随机场的离散方法与有限元法的离散方法类似,按照离散后随机变量与离散单元随机场的关系,随机场的离散方法主要有中心点法、局部平均法、形函数插值法、加权积分法和正交展开法等。局部平均法是利用随机场单元中的平均值、均值方差来代表随机场单元的性质。局部平均法具有使用方便、计算精度高,且对随机场相关结构不敏感等特点,而成为一种最为常用的随机场处理方法。自局部平均法理论提出之后,刘宁[6]提出三维可分向量随机场局部平均法,能够较好地离散三维随机场。李永见[7]提出普遍适用于岩土工程土性参数随机场的三维随机场离散模型,既考虑了随机场的三维特性,又便于实现。本文基于以上两种方法对渗透系数空间三维随机场进行局部平均离散,进而得到一系列渗透系数随机变量。渗透系数的空间分布通常为对数正态分布[8,9],本文渗透系数随机场取为对数正态分布的平稳随机场,即空间内各点渗透系数均值相同,两点之间渗透参数的相关性只与两点之间的距离相关。

设α(x,y,z)(x,y,z∈V)为三维连续的宽平稳随机场,根据宽平稳随机场的定义有:

(1)“点”均值函数为常数:

(2)“点”方差函数为常数:

(3)自协方差函数只与相关步距有关:

则单元内的局部平均为:

在有限元计算中,式(1)可进一步表示为:

式中:Ve——单元e的体积;

αe(ξk,ηk,ζk)——随机场α(x,y,z)在第e单元第k高斯点处的值;

ξk、ηk、ζk——分别为第k个高斯点的局部坐标;

Wk——高斯点加权系数;

|J|k——Jacobi矩阵行列式值;

ng——单元内高斯点数目。

由式(2)可得单元内局部平均随机变量的均值和协方差为:

其中:

式中:——第e单元第k高斯点和第e'单元第l高斯点的互协方差;

——分别为两点在x、y、z方向的距离;

ρw——相应的相关函数值。

应用三维可分离随机场法,相关函数可分离为:

则任意两个局部平均随机场单元e、e'协方差的高斯积分表达式为:

式中:i=j=k=p=q=r=8为积分点个数;

(ξi,ηj,ζk)、)——分别为各积分点的局部坐标;

Wi、Wj、Wk、Wp、Wq、Wr——分别为加权系数。

应用此原理在程序计算时,首先计算各积分点处的形函数及形函数对局部坐标的偏导,然后计算形函数与单元各节点坐标之积,代入相关函数表达式,计算相关函数在各积分点处的值,再依次求出各单元在各积分点处的Jacobi矩阵行列式值,最后通过六重求和计算出随机场任意两单元之间的协方差,并组合成协方差矩阵[C]n×n。模型中各单元为任意六面体单元,协方差矩阵[C]n×n由式(9)计算,由于积分点个数i=j=k=p=q=r=8,计算量较大,计算速度较慢。在利用公式(9)计算时,要预先给定标准自相关函数的表达式。Sudicky(1986)通过试验研究得出渗透系数的自相关函数可以用指数函数形式描述[10]:

式中:δ1、δ2——分别是水平向和垂直向上的相关距离。

3 三维渗流摄动随机有限元法的原理

三维稳定渗流的有限元列式为:

式中K为总体渗透矩阵,依赖于介质的渗透特性和节点几何特性。如果介质的渗透性能有随机扰动,则此随机扰动可作为空间的随机场或随机过程,这样渗透矩阵K就是随机的。另外节点等效流量向量F也可能是随机扰动的。由控制方程可看出,节点水头列阵H(即系统响应)必定具有随机扰动性。将渗透张量看作各向异性的随机场。假设αi为渗透系数ki在均值点处的微小摄动量,则,运用中心二阶摄动方法可得到如下递推方程组:

式中:K0、H0、F0——分别表示渗透系数、节点水头、节点水头外荷载在均值点的值,下标i、j表示对αi、αj求导。

依次求解上述递推方程组可得水头的一阶和二阶变异量,由此得到节点水头均值为:

节点水头协方差为:

利用以上计算结果,进一步可得水力梯度的均值和协方差如下[11]:

当随机变量总数N较大时,计算时间较长,为此本文程序引入特征正交化技术减小计算量。由特征正交化可知,N维相关正态分布随机向量α={α1,α2,…,αN}可以用独立正态分布的随机向量β={β1,β2,…,βn}来表示,则:

式中:Φ——协方差矩阵的特征向量矩阵。

对于两个随机变量βi、βj,有:

式中:λi为协方差矩阵的特征值。

式(13)、(14)、(15)可表达为:

其中:

由此可得水头均值和协方差为:

本文应用局部平均法对三维渗透系数随机场进行离散,同时应用特征正交化对摄动有限元进行简化,根据上述理论,笔者编制三维渗流摄动随机有限元程序,具体程序流程图如图1所示。

4 三维渗流摄动随机有限元计算程序数值论证

为验证本文所编制的三维渗流摄动随机有限元计算程序的正确性,应用一个典型算例,分别采用三维渗流摄动随机有限元法和Monte Carlo随机有限元法计算,并进行比较。数值验证模型断面如图2所示。边界ab、cd上水头已知,分别为h1=10m、h2=2m,且两个边界上水头为确定值,不考虑水位的波动,在计算时水头的变异系数为0。模型土体为均质土,渗透系数的均值为1.16×10-5m/s,标准差为2.31×10-6m/s,相关距离为3m,相关函数取为三维高斯型相关函数。对渗流场进行有限元网格和随机场网格剖分时,采用相同的网格,网格为1.0m×1.0m×1.0m的六面体网格。

通过两种方法计算的渗流场水头均值的等势线如图3所示。从图3可看出,两种方法计算的水头均值基本相等,相对误差在1%以内。两种方法计算的水头标准差分布如图4所示。由图4可看出,两种方法计算的水头标准差也基本相等,相对误差在5%以内。

5结语

本文将渗透系数看作空间随机场,应用局部平均法对渗透系数随机场进行离散,编制三维渗流摄动随机有限元程序,计算随机渗流场。通过与蒙特卡罗法计算的水头均值和标准差比较可知,本文编制的三维渗流摄动随机有限元程序计算的随机渗流场与Monte Carlo随机有限元法计算的渗流场相近,由此证明了本文所编制的三维渗流摄动随机有限元程序的正确性和可行性。

摘要：在非确定性渗流研究中,本文考虑渗透系数的空间随机性,基于摄动随机有限元法,推导相应的渗流场随机分析中渗流响应量(水头)的随机响应公式,进而实现三维稳定渗流场随机有限元分析。通过一个算例,将三维摄动随机有限元法计算的随机渗流场与蒙特卡罗法计算的渗流场进行对比,验证自编三维随机渗流摄动有限元程序的正确性和可行性。

关键词：摄动随机有限元,蒙特卡罗法,渗流计算

参考文献

[1]Griffiths D V.Three-dimensional seepage through spatially random soil[J].Journal of Geotechnical and Geoenvironmental Engineering,1997,123(2):153-160.

[2]Handa K,Anderson K.Application of finite element methods in statistical analysis of structures[C].Proc3rd Int Conf on Struct Safety and Reliability Trondheim,Norway,1981,409-417.

[3]Hisada T,Nakagiri S.Stochastic finite element analysis of uncertain system[C].Proc 4th Int Conf on FEM,Auatrulia,1982,133-138.

[4]Erik Vanmarcke.Random Fields;Analysis and Synthesis[M].Cambridge:The MIT Press,1983.

[5]陈虬,刘先斌.随机有限元及其工程应用[M].西安:西安交通大学出版社,1993.83-95.

[6]刘宁.三维可分向量随机场局部平均的三维随机有限元及可靠度计算[J].水利学报,1995,(6):75-82.

[7]李永见,王月香,茜平一.土性参数随机场的离散研究[J].四川建筑科学研究,2002,28(2):37-39.

[8]施小清,吴吉春,袁永生.渗透系数空间变异性研究[J].水科学进展,2005,(2):210-215.

[9]束龙仓,李伟.北塘水库库底底层渗透系数的随机分析[J].吉林大学学报(地球科学版),2007,37(2):216-220.

[10]Sudicky E A.A natural gradient experiment on dolute transport in a sand aquifer;Spatial variability of hydraulic conductivity and its role in the dispersion process[J].Water Resour.Res,1986,22(13):2069-2083.

三维有限元法 篇11

水库工程位于广西的灵奇河源头坤屯河上, 是一座以发电为主、兼有防洪等综合效益的枢纽工程。坝址以上集雨面积110.7 km2, 总库容3747.26万m3, 多年流量3.00 m3/s, 多年平均径流量9475.6万m3。工程由挡水建筑物、引水式电站、引水渠道及渠系建筑物组成, 挡水建筑物混凝土双曲拱坝。坝高65.5 m, 坝顶上游面半径126.6 m, 坝顶上游面弧长245.1 m, 坝顶弧长208.5 m, 坝顶中心角110.9º, 底拱中心角59º, 坝底厚15.5 m, 顶厚4.0 m, 非溢流坝顶高程586.50 m, 溢流堰顶高程为580.0 m。采用表孔高低坎差动式挑流消能, 溢流前缘净宽42.5 m。

水库于1996年建成拱坝, 2008年相关的设计研究机构对水库大坝进行了安全评价。因坝后工作栈道未完建, 溢洪道、左侧坝肩山体、进库及坝顶交通、坝后消力池、监测系统等尚未完善配套, 枢纽工程也未进行完工验收, 不具备向沿途乡镇乃至县城供水的功能, 大坝为二类坝。因此需对水库进行续建和配套, 尽早发挥枢纽的综合效益,

二、水库扩容方案

水库因不具备加高的地形地质条件, 扩容主要针对溢流坝段进行。将敞泄式溢流堰改建为控制泄流型式, 提高正常蓄水位3 m、扩大有效库容, 扩容前后特征水位表1所示。具体方案为将上游弧形悬挑部分从顶部打去0.7 m, 即溢流堰顶高程降至579.3 m, 顶部平台两边做成弧形衔接以利过水, 前端悬挑部从下部补强;溢流堰顶加设5孔露顶式液压平面闸门。大坝正常蓄水位抬高后, 拱坝的应力发生变化。为了解扩容后拱坝的应力状况, 本文将采用大型有限元软件ANSYS对其位移和应力进行计算分析, 为设计提供依据。

三、计算模型建立

1. 模型尺寸

采用三维实体模型, 整体坐标为OXYZ。采用右手坐标系, X轴向左岸为正, Y向河道下游为正, Z向上为正。模型底部为固定边界, 上下游端侧面为Y方向约束, 左右端面为X方向约束。为保证结果的精确性和后处理的方便性, 在不能保证所有网格为六面体时的情况下, 优先保证坝体为六面体网格, 共有单元数23 310个, 节点数103 986个。根据拱坝的实际情况, 将拱坝分为9层, 上下游坝面曲线均采用二次抛物线。考虑新增水闸的制作, 形成一个长×高为44.5×3.5 m的缺口, 模型网格图如图1、图2所示。

2. 材料物理力学参数

坝体由2种混凝土组成, 基础由6种岩体材料组成, 均为各项同性线弹性材料, 其材料参数如表2所示。

3. 计算工况

计算工况及荷载主要由工程前2种工况和工程后4种工况。工况1:工程前, 原正常蓄水位+温降工况, 荷载为坝体自重+水压力+淤沙压力+温降温度荷载;工况2:工程前, 原正常蓄水位+温升工况, 荷载为坝体自重+水压力+淤沙压力+温升温度荷载;工况3:工程后, 正常蓄水位+温降工况, 荷载为坝体自重+水压力+淤沙压力+温降温度荷载;工况4:工程后, 正常蓄水位+温升工况, 荷载为坝体自重+水压力+淤沙压力+温升温度荷载;工况5:设计洪水位+温升工况, 荷载坝体自重+上游水压力+下游水压力+淤沙压力+温升温度荷载;工况6:工程后, 校核洪水位+温升工况, 荷载为坝体自重+上游校核洪水位时水压力+下游校核洪水位时水压力+淤沙压力+温升温度荷载。其中, 封拱温度为18℃, 采用Tm1=Tm0, Td1=Td0。坝体温度场可以按与应力分析相同的有限元网格模型进行计算得到。

四、计算结果分析

1. 坝体位移分析

由表3可知, 各工况下, 坝体位移的最大值均位于拱冠处。工程增效扩容前, 坝体最大位移正常蓄水位+温降工况控制, 位移为34.16 mm, 工程后最大位移为40.65 mm, 位移增加了6.49 mm;正常蓄水位+温升工况时, 工程增容扩效后, 位移量为15.58mm, 增加了3.69 mm;工程扩容后, 设计洪水位+温升及校核洪水位+温升工况下, 最大位移量分别为17.42 mm、20.26 mm。

2. 坝体应力分析

坝体主拉和主压有限元等效应力如表4所示。

正常蓄水位+温降工况下, 温降不仅使上游面边缘产生较大的拉应力, 也使下游面中部产生较大的拉应力。扩容前, 上游面最大拉应力出现在坝体右端高程525 m处, 应力值为0.92 Mpa, 最大压应力出现在拱冠高程555.8 m处, 应力值为-2.69 Mpa;下游面最大拉应力出现在坝体左端高程525 m处, 应力值为1.09 Mpa, 最大压应力出现在拱冠高程525 m处, 应力值为-1.59 MPa。扩容后, 上下游面的拉压应力值均有所增加, 上游面最大拉应力为1.4 Mpa, 增加了0.48 MPa;最大压应力值为-3.13, 变化值为0.44 MPa;下游面最大拉应力值为1.47 Mpa, 增加了0.39 Mpa, 最大压应力值为-2.35 Mpa, 变化了0.76 Mpa。

正常蓄水位+温升工况下, 温升不仅使下游面中部较大的拉应力消失, 也使上游面边缘拉应力逐渐减少甚至转变为受压。扩容前, 上游最大拉应力位于坝体左端处高程576 m, 应力值为1.16 Mpa, 最大压应力值位于拱冠高程566 m处, 应力值为-1.51 Mpa;下游面最大拉应力位于坝体左端高程525 m, 应力值为1.17 Mpa, 最大压应力位于拱冠高程545.5 m, 应力值为-2.58 Mpa。扩容后, 上游面最大压应力值为1.37 MPa, 增加了0.21 MPa, 最大压应力值为-1.98 MPa, 变化值为0.47 MPa;下游面最大拉应力值为1.30 MPa, 增加了0.13 MPa, 最大压应力值为-2.71 MPa, 变化值为0.13 MPa。

扩容后, 设计洪水位+温升工况下, 上游面最大拉压力为1.45 MPa, 最大压应力为-2.2 MPa, 下游面最大拉应力为1.37 MPa, 最大压应力为-2.78 MPa;校核洪水位+温升工况下, 上游面最大拉压力为1.55 MPa, 最大压应力为-2.49 MPa, 下游面最大拉应力为1.46 MPa, 最大压应力为-2.87 MPa。

由各种工况计算成果可推断, 温降是拉应力控制工况, 温升是压应力控制工况。最大等效拉应力为1.47 (工况3) <[1.50 MPa], 最大压应力为-2.866 MPa (工况6) , 其绝对值小于非地震情况的特殊组合的容许压应力[5.5 MPa]。

五、结语

本文采用有限元法对水库拱坝扩容前后的坝体位移及应力状况进行了研究, 得到以下几点结论。

(1) 坝体最大位移由工程后的正常蓄水位+温降工况控件, 最大位移处为坝顶拱冠处:40.65 mm, 工程后正常蓄水位+温降和温升工况的坝体最大位移分别增加了6.49 mm和3.69 mm, 位移绝对值和增加值都较小。

(2) 温度降低使上游面边缘产生较大的拉应力, 下游面中部产生较大的拉应力;温度升高使下游面中部较大的拉应力消失, 上游面边缘拉应力逐渐减少甚至转变为受压。

(3) 水库扩容后, 各工况下, 坝体上下游面的主拉及主要应力极值均有所变大。温降工况下, 最大拉应力增加了0.48 MPa, 最大压应力增加了0.76 MPa;温升工况下, 最大拉应力增加了0.21 MPa, 最大压应力增加了0.47 MPa。

(4) 由各种工况计算成果可推断, 温降是拉应力控制工况, 温升是压应力控制工况。最大等效拉应力为1.47<[1.50 MPa], 最大压应力为-2.866 MPa。其绝对值小于非地震情况的特殊组合的容许压应力[5.5 MPa], 满足规范要求。

参考文献

[1]孙大为, 蔡海勇等.长岭水库拱坝加固前后三维有限元分析[J].人民黄河, 2010 (10) .

[2]SL282-2003, 混凝土拱坝设计规范[S].北京:中国水利水电出版社, 2003.

[3]赵玮.拱坝三维有限元等效应力计算研究[D], 2007.

[4]李同春, 陈回芳等.网格尺寸对拱坝等效应力分析的影响[J].水利学报, 2004 (09) .