仿真数值计算（精选9篇）
仿真数值计算 篇1
1工程基本情况[1]
富川水电站位于云南省维西傈僳族自治县境内澜沧江右岸一级支流其普河下游,是一座径流式水电站,坝址集水面积87.75 km2,多年平均径流量3.81 m3/s,总装机40 MW,年发电量22 585万kWh。电站拦河坝为坝高27.5 m的混凝土埋块石重力拱坝,正常蓄水位2 256.0 m,由一条长6 050 m的隧洞和1 500 m长的压力钢管引水进入厂房,钢管外径从1.6～1.4～1.3 m逐级变化,最大发电流量8.08 m3/s,发电净水头623 m。
2岔管设计布置参数
富川水电站岔管靠山边布置在发电厂房后,卜形月牙肋岔管,分岔角55°,为空间平面岔管,其平面与水平面夹角为42°,岔管中心高程1 635.7 m。岔管公切球半径800 mm,壁厚38 mm,肋板厚78 mm;岔管前的主管内径1 228 mm,壁厚36 mm;岔管后的支管内径848 mm,壁厚30 mm,岔管体型及结构布置见图1、图2。设计采用的岔管材料为WDB620高强钢,弹性模量2.1×105 N/mm2,泊松比0.3,密度7 850 kg/m3。
3结构数值模型仿真计算
3.1岔管结构有限元仿真计算
为了给富川水电站岔管结构设计提供可靠依据,研究人员对岔管结构用有限元方法进行仿真计算,主要计算该岔管在其所承担荷载工况(包括设计静水压力工况、运行工况和试验工况)下的应力、位移,分析岔管应力集中部位的应力分布情况,设计人员可以据此校核结构应力是否满足材料强度,并评价岔管在不同荷载工况下的应力安全储备情况。
3.2荷载组合情况
该岔管静力有限元计算各工况取值:
(1)正常运行工况。正常蓄水位下静水压力623×1.1=685.3 m,约6.72 MPa。
(2)特殊运行工况。校核洪水位下静水压力627.6×1.1=690.36 m,约6.77 MPa。
(3)试验工况。试验压力水头为623×1.1×1.25=856.63 m,约8.40 MPa。
在进行有限元计算时,上述静水压力水头均以压力强度P=γ0h施加,h为计算水头。
3.3计算假定[2]
富川水电站岔管有限元模型见图3。为使建立的力学模型在符合荷载及边界条件下简化计算,便于离散方案的实现,作如下基本假定:①岔管结构假定单独承担内水压力,不计岩体或混凝土抗力,即按明岔管进行有限元计算。②假定岔管结构符合线弹性基本假定,结构按线弹性有限元计算。③计算荷载按静水压力模拟,对于水击(锤)压力按等效静水压力计算。
3.4计算方法[2]
富川水电站的设计研究人员采用ANSYS有限元计算软件,对于月牙肋钢岔管这样的复杂空间结构进行有限元划分时,遵循以下原则:① 计算模型应尽量与实际运行的岔管在结构特征方面取得一致;② 所选用的离散单元应尽量模拟岔管受力和变形特征;③ 荷载和边界条件应接近岔管的实际工况。
4研究计算成果分析[3]
4.1正常运行工况结果分析
4.1.1位移分析
(1)整体管壳位移分析。
从有限元计算结果来看,在该工况荷载作用下,由于岔管为“卜”形结构,岔管管壳位移成不对称分布,呈“鼓包”形态(见图4、图5);最大总位移值(USUM值)为1.1 mm,位于左支锥管外壁(图5),位移云图中的总位移值依次递减;整体管壳总位移沿水流方向平面上下基本对称(图5);支锥管上下表面表现为外凸形态(±Y方向),最大外凸位移值0.9 mm;左支管侧壁外侧(截面1 800左右位置)呈内凹形态(X方向),右支锥管侧壁外侧(截面00左右位置)呈外凸形态(X方向),最大位移值0.7 mm。
(2)月牙肋板位移分析。
从计算结果来看,由于结构的非对称性,月牙肋板位移特征表现为不仅Y方向有变形,侧向(X方向)也有变形,最大Y方向位移为0.38 mm,位于月牙肋板下肢截面最小处(图6);最大侧向(X方向)位移为0.30 mm,位于上、下肢截面最小处;月牙肋板总位移为0.58 mm,位于上、下肢截面最小处(图7)。从月牙肋板的变形特征来看,肋板呈现双向变形特点,由此可得出其应力特征:在肋板截面最大处将出现最大应力;由于肋板的侧向变形以及支锥管的外凸变形作用,肋板与左、右支锥管连接处(岔裆部位)将出现最大主应力,且均为拉应力。
4.1.2应力分析
为减小边界约束对岔管结构应力的影响,将主、支管沿水流方向延伸1倍管径后施加边界约束,在约束边界附近的管壁属于边界条件影响区域,因此进行应力分析时,取出远离边界条件影响区域的岔管结构进行分析,即锥管壁厚度为34～38 mm的岔管结构部分及月牙肋板。根据规范规定[4],钢岔管各计算点应力应按第4强度理论进行计算校核。
(1)岔管结构管壳应力分析。
①整体膜应力。岔馆整体膜应力区域主要是管壁远离加强构件部分的区域,从计算结果来看,应力分布特征呈与水流方向平面对称分布,即上、下表面分布基本相同。膜应力区最大膜应力出现在左支锥管处,其最大值约160 MPa,沿四周逐渐减小,管壁膜应力值在150 MPa以上的区域很小,见图8。
②局部膜应力区应力。局部膜应力区主要出现在管壁突变以及月牙肋板近旁管壁处。岔管结构体型突变处主要是指以下部位:主管与主锥管连接处、相关线交汇处、锥管母线转折处、支直管与支锥管母线转折处等部位。从计算的结果来看,结构体型突变处的应力较其他部位大,如母线转折部位(主要是左支锥管)主、支锥管母线交汇处(E点内壁)较其他体型突变的部位高,最大主应力及Von Mises应力为214 MPa,其他部位的应力值在190 MPa以内,局部膜应力分布如图9所示。
③局部膜应力+弯曲正应力分析。在岔管月牙肋与管壁连接处及其局部范围内(岔裆部位)出现局部膜应力和弯曲正应力月牙肋与管壁连接处靠近月牙肋肢尖附近应力较大,其附近应力最大值为204 MPa,最大Von Mises应力出现在月牙肋上肢尖与管壁连接处,其值为236 MPa。
(2)月牙肋板应力分析。
从计算结果来看,由于岔管整体结构的不对称,月牙肋板呈双向变形特征,其应力分布具有沿厚度方向分布不对称的特点,但沿上、下两肢方向应力分布较为对称。最大Von Mises应力出现在月牙肋板截面最大处,其值为174 MPa,Von Mises应力值沿月牙肋两肢方向逐渐递减。
4.2特殊运行工况结果分析
4.2.1位移分析
从有限元计算结果来看,在特殊运行工况下,岔管管壳和月牙肋板位移的分布形态均与正常运行工况基本相似。
(1)整体管壳位移。
最大总位移值为1.2 mm,位于左支锥管外壁;支锥管上下表面最大外凸位移值0.86 mm;左支管侧壁外侧(截面1 800左右位置)呈内凹形态(X方向),右支锥管侧壁外侧(截面00左右位置)呈外凸形态(X方向),最大位移值0.66 mm。
(2)月牙肋板位移。
最大Y方向位移为0.38 mm,最大侧向(X方向)位移为0.30 mm,月牙肋板总位移为0.58 mm。
4.2.2应力分析
从计算结果来看,在特殊运行工况下,岔管结构管壳和月牙肋板应力分布形态与正常运行工况基本相似。
(1)岔管结构管壳应力分析。
①整体膜应力。膜应力区最大膜应力出现在左支锥管处,膜应力最大值为161 MPa,沿四周逐渐减小;应力值在150 MPa以上的区域范围较小。
②局部膜应力区应力。计算结果表明,结构体型突变处的应力较其他部位大,如母线转折部位(主要是左支锥管)主、支锥管母线交汇处(E点内壁)较其他体型突变的部位高,最大主应力及Von Mises应力为215 MPa,其他部位的应力值在180 MPa以内。
③局部膜应力+弯曲正应力。月牙肋与管壁连接处靠近月牙肋肢尖附近应力较大,其附近应力最大值为211 MPa,最大Von Mises应力出现在月牙肋上肢尖与管壁连接处,其值为238 MPa;大于150 MPa以上的应力区域较小,局部应力由于应力集中只在局部范围内出现。
(2)月牙肋板应力分析。
计算结果得出最大Von Mises应力出现在月牙肋板截面最大处,其值为176 MPa,Von Mises应力值沿月牙肋两肢方向逐渐递减,Von Mises应力值大于150 MPa应力分布范围很小。
5结语
由三维有限元数值仿真计算结果,设计技术人员对岔管进行了强度复核,在岔管采用WDB620高强钢时,其在各个工况下的应力均满足规范规定,可以满足运行要求。总结富川水电站岔管设计,如采用过去近似计算方法,对体型较为复杂的不对称钢岔管分析时往往非常麻烦,且计算成果不甚满意,而采用有限元技术进行数值仿真计算则较为直观,对钢岔管的设计而言不失为一种全新途径。
摘要：富川水电站引水岔管为“卜”形月牙肋,属高水头、大HD值的岔管结构。设计、研究人员合作对岔管结构进行三维有限元数值仿真计算,研究分析在各种荷载组合下岔管结构管壳应力、月牙肋板应力,为岔管结构安全设计提供依据。
关键词：富川水电站,岔管设计,仿真
参考文献
[1]上饶市水利电力勘测设计院.云南省维西傈僳族自治县富川水电站可行性研究报告[R].2010.
[2]南昌大学,上饶市水利电力勘测设计院.富川水电站岔管静力三维有限元仿真模拟分析计算研究报告[R].2011.
[3]谢冠峰,李火坤.高水头卜型月牙肋岔管结构三维应力分析[J].南昌大学学报(工科版),2009,(1):90-95.
[4]SL281-2003,水电站压力钢管设计规范[S].2003.
仿真数值计算 篇2
航空母舰尾流数值仿真研究
为研究舰载飞机进舰着舰过程中航空母舰尾流的影响,选取了一个尾流模型在微机的MATLAB环境中模拟其速度的空间分布.根据成因将紊流分为4个组成部分:自由紊流分量的.功率谱用成形滤波器法模拟,稳态分量通过线性插值得到,周期性分量可直接计算得到,而航母诱导的随机紊流分量则通过MATLAB的SIMULINK工具包模拟.模拟表明周期性分量占主导地位.此模型需要实验数据和船上实测数据进一步确认.
作 者：彭兢 金长江 PENG Jing JIN Chang-jiang 作者单位：北京航空航天大学,飞行器设计与应用力学系刊 名：北京航空航天大学学报 ISTIC EI PKU英文刊名：JOURNAL OF BEIJING UNIVERSITY OF AERONAUTICS AND ASTRONAUTICS年，卷(期)：26(3)分类号：V211.78关键词：航空母舰 舰载飞机 尾流 紊流 功率谱
仿真数值计算 篇3
气体放电产生的低温等离子体在材料表面改性方面有许多重要的应用,并显示出很好的经济效益,尤其是环保效益[1,2]。对于等离子体表面处理在工业上的应用来说,最佳的选择是功率适中(几百mW/cm3)、放电均匀的大气压下辉光放电(APGD)产生的等离子体。利用氦气比较容易实现大气压辉光放电。
从上个世纪80年代末开始,APGD成为人们的研究热点[4] ,通常采用介质阻挡的平行板电极结构。由于人们急需摸清大气压辉光放电的放电机理,而纯粹的实验不能完全提供必要的内部过程,因此对大气压辉光放电的数值模拟计算成为当前这个领域内的研究热点。法国的Massines小组较早地对氦气和氮气的APGD进行了实验研究和数值模拟[5],他们认为:得到大气压下均匀放电的关键是在较低电场下缓慢发展大量的电子雪崩。明尼苏达大学研究小组随后也进行了氦气APGD的实验研究和数值模拟[6],他们在计算中计入了亚稳态与彭宁效应的影响。大连理工大学的王德真也进行了纯氦气APGD的数值模拟[7],并取得了一定的成果。然而,以往的数值模拟计算中,由于人们对APGD放电机理认识不足、参数选取混乱等诸多原因,导致常常出现计算结果与实验结果不能很好地对应,比如计算得到的放电电压远高于实际中放电所需要的。这就大大降低了利用仿真的计算结果来代替实际放电过程中不易测量的内部参数的可信度,这也是APGD放电机理研究所遇到的困难之一。
本文细致地选择了放电计算中所使用的各个参数,力图确保它们最接近实际值而不是理论值,这使得最终的计算结果与同步的实验结果有很好的吻合度。这样的计算结果将对实验具有就较高的实用价值,对进一步深化了解大气压氦气辉光放电的机理具有十分重要的意义。它揭示了放电过程中许多不易测量出的参量随时间、距离的演化,并且本文中已经使用这些计算出的结果对实验中观察到的现象做出了理论解释。
2 理论模型
2.1 电路结构
由于在介质阻挡放电中,平板电极的直径远大于气体间隙的距离,因此可采用一维流体模型并假设:间隙距离为d,上介质下表面为x=0,是第一个离散格点N=1;下介质上表面为x=d,是最后一个离散格点N=L。电源电压为Vapp(t)=Vmsin(2πft),下电极接地,外加电场参考方向为x方向。
2.2 基本方程
介质阻挡辉光放电的一维流体模型由离子、电子连续方程和泊松方程构成。对于大气压氦气辉光放电,模型中还作了以下的假设:
He是非电负性气体,不考虑电子附着到原子(分子)上而形成负离子;只考虑天然本底辐射导致的电离G0以及He和电子的直接碰撞电离所产生的e和He+,忽略离子转换而生成He分子离子[8]以及重分子离子;忽略气体中的杂质及其电离,尤其是亚稳态He(23S)和N2的彭宁电离;只考虑e和He+复合而导致的带电粒子损失;由于是一维模型,忽略了径向扩散导致的带电粒子损失。
在这些假设和简化下,有如下的方程组:
主方程:
undefined
辅助方程:
undefined
其中,q代表元电荷;p,n分别代表放电空间中的离子密度和电子密度;ψ为放电空间中的电位(假定接地点为0电位),E为电场强度;G代表电子离子产生项,L代表电子离子复合项;Jp,Jn分别代表离子和电子在气体间隙中运动的运流电流密度;μp,μn分别代表离子迁移率和电子迁移率,数据都来自文献[9];Dp,Dn分别代表离子扩散系数和电子扩散系数,数据分别来自文献[11,12];α代表电子碰撞电离系数,数据来自文献[11];r代表电子-离子复合率,数据来自文献[10]。
实验中绝缘介质采用石英玻璃,其相对介电常数为εr=3.8,单片介质厚度为0.1cm,放电间隙d=0.5cm,电极面积S=19.63cm2,外加电源电压幅值为1300V,频率为30kHz。
2.3 边界条件与初始条件
泊松方程的边界条件由介质与气体分界面的电磁场边界条件给出:
undefined
其中Eg、Ed分别为气体和介质中的电场强度,εd为介质的介电常数,σ为介质表面积累的面电荷密度。
对于连续方程,在边界处考虑二次电子发射,有:
undefined
其中,Jc为总运流电流密度,γ为二次电子发射系数,本文中γ取0.1。
放电空间中的初始电子主要来自于天然本底辐射。由中国辐射防护研究院对充满氦气的封闭钢罐所作的测量,有:G0=4.95对(离子、电子)/s/cm3。由于放电前气体间隙中带电粒子的产生和消失处于平衡状态,因此有G0=rneni=rnundefined,求得初始电子密度为:n0=1.7×102/cm3。
3 模拟结果与分析
本文采用Scharferter和Gummel提出的电流离散形式和有限差分方法[13]对模型中的非线性方程组联立求解。
图1给出了放电稳定后两个周期的放电仿真图形。每半个电压周波就有一次放电,并且具有很稳定的周期性。大约每半个周期开始4微秒后放电电流达到峰值,为30毫安。此电流为全电流,即运流电流与位移电流之和。
图1给出了外加电压、介质压降、气体压降、以及放电电流这些电路外特性参量随时间的演化过程。外加电压从0V开始上升,气体电压也随着迅速上升。随着气体电场增强超过一定值,气体内电离作用迅速增强,产生大量带电粒子,导致电流剧增,于是气体状态由绝缘变为导通。由于存在介质阻挡,介质等效电容由于电路中总电流的迅速增长而迅速积累与外加电压同向的电压,这一机制导致气体上的压降从1200V左右迅速跌落至500V以下,气体中的电场也因此迅速跌落,强烈的电离过程被终止。随后气体中的带电粒子在迁移、扩散和复合机制的共同作用下迅速消失,气体恢复绝缘,放电电流也迅速下降。而介质上的电压仍保持高位,当外加电压进入下半个周期反向后,气体上的电压等于外加电压与介质电压的叠加,于是气体电压比外加电压提前达到1200V左右的击穿电压值,发生第二次击穿,并很快由于气体电压减小而恢复绝缘,然而此时电压还未达到峰值,于是气体电压在下落后还会随外加电压有小幅上升,大约升至500V左右。整个放电过程说明了介质阻挡机制的引入对实现稳定的辉光放电的重要性。
图2为与仿真计算采用相同参数的实验波形。通过对比可以看出,仿真结果与真正的实验波形吻合得很好:外加电压都在1300V左右,放电电压都在1200V左右,放电之后气体电压峰值也都在500V左右,放电电流都在30mA左右。这说明仿真结果是对实验波形很好的近似,这样的仿真结果为实验提供了有较高的可信度的内部参数。
图3给出了气体中的电场在半个周期内随时间的演化过程。图中所对应为正半周放电,因此x=0cm处为阳极,x=0.5cm处为阴极。可以看到电场在放电时呈现很明显的阴极位降区,阴极处的电场峰值达到接近10kV/cm。在放电结束后由于气体电压迅速下降使得电场迅速恢复到各处比较均匀的状态,并且随气体电压再次小幅上升而上升。
图4为在与仿真计算采用相同参数的实验中用ICCD相机在放电不同时刻所拍摄的照片:上半部分为拍摄时刻与电流大小的对应关系,下半部分为对应时刻的短曝光时间的放电照片。图5为在相同时刻的放电内部参数的仿真计算结果,六幅图按时间发展排序为从左至右,从上至下。
在图5第一、二张图这两个阶段气体空间中电场分布还较为均匀,电子在整个间隙中拥有的能量差别不大,因此激发态粒子密度及其辐射光强正比于电子密度,而计算得到的电子密度正是在阳极处最大,并且随时间发展阳极的电子密度进一步加大。这正好解释了图4(a)中看到的放电区在阳极区有微弱的发光,并且(b)中阳极发光得到增强。这意味着此时的放电仍为汤森放电过程。图5第三张图中,间隙中电场已经发生畸变,整个间隙中离子密度和电子密度都很高,靠近阴极部分已经进入1010/cm3数量级。这对应着图4(c)放电空间全亮并且靠近阴极部分亮度更大。这标志着随着电压增大,放电正在由汤森放电转化为辉光放电。图5第四张图具有明显的辉光放电内部参数特征,图4(d)也具有明显的辉光放电外部发光特征,说明放电已经进入辉光放电阶段。图5第五、六张分别展示了放电电流最大时刻和开始熄灭时刻的内部参数变化过程,阴极位降区厚度不断地减小,电子离子密度也很快减少,对应着图4(e)、(f)中负辉光区变薄与整个放电空间发光减弱。
图6给出了放电电流最大时刻气体中的离子、电子密度和电场随距离的分布,大气压氦气辉光放电具有十分类似低气压辉光放电的特点:即存在阴极位降区、负辉光区、法拉第暗区和正柱区。在本文所采用的放电条件下,阴极位降区的宽度大约0.5mm,而对应的实验所测定的负辉光区宽度也为0.5mm;正柱区大约2mm左右,其中电子、离子密度相等,大约在1010/cm3左右。
4 结论
基于燃面推移的内流场数值仿真 篇4
基于燃面推移的内流场数值仿真
基于处理移动边界的动网格技术,采用隐式的、非定常的耦合解法,选用k-ε双方程的Realizable湍流模型和壁面函数法,对某小型固体火箭发动机粘性的、可压缩的、高雷诺数的燃烧室-喷管湍流流场进行了一体化的数值计算,得到其在移动边界条件下瞬态的.流场分布,并分析了其完整的内弹道各参数随着装药燃面推移而连续变化的规律,通过数据拟合出的三条内弹道特性参数曲线与试验数据符合较好,说明文中所采用的内弹道数值仿真计算方法可靠.
作 者：谢丽宽 马拯 俞红博 于冀国 XIE Li-kuan MA Zheng YU Hong-bo YU Ji-guo 作者单位：海军驻西安地区军事代表室,西安,710025刊 名：弹箭与制导学报 PKU英文刊名：JOURNAL OF PROJECTILES, ROCKETS, MISSILES AND GUIDANCE年，卷(期)：27(3)分类号：V435关键词：固体火箭发动机 数值模拟 燃面
盾构的数值仿真分析 篇5
1.1 ANSYS简介
ANSYS是目前广泛应用于机械工程领域的数值分析软件,通过在计算机上模拟现场的各种条件,通过节点间的连接模拟构件的实际受力情况,通过加载受力,在计算机上提前预测各种不利效应,从而判断该结构是否满足安全、使用要求。
1.2 ANSYS软件主要功能
ANSYS软件能模拟大部分物理力学作用,可广泛的用于地下工程、航空航天、机械摩擦和制造、物理碰撞、桥梁施工、温度效应、大型水利工程、等诸多专业学科。能对结构进行自适应单元划分,从而进行非线性分析,通过对参数的调整,确定各影响因素的敏感性,对各种设计方案进行优化,从而在实际应用中,通过各种物理化为作用,对敏感度高的影响因素重点改进,对耗费经济少,功效较高的影响因素重点优化,确定最优设计方案。
在盾构施工中,由于土层的受力比较复杂,根据经验难以断定地层的损失量,采用过于保守的设计会加大工程造价,通过ANSYS模拟盾构隧道施工全过程,可以有效的控制工程造价,保证地下工程施工安全。
2 盾构隧道分析
采用四结点单元划分网格如图2所示:
首先我们研究盾构开挖时所产生的沉降变形,由于衬砌单元具有很强的刚度,开挖过程中并没有衬砌单元的存在,于是先采取单元生死技能杀死衬砌单元便可以得到,要激活“单元死”的效果,ANSYS程序并不是将“杀死”的单元从模型中删除,而是将其刚度矩阵乘以一个很小的因子[ESTIF]。从而不对载荷向量生效。
从云图4可以看出,地表沉降在隧道盾构正上方处沉降量最大,离正上方水平距离越远,沉降越来越小,成正态分布,下面采用peck理论公式计算进行对比。
Peck公式为:Sx=Smaxexp(-x2/(2i2))
式中:Sx——x处的沉降量。
x——离中线的水平距离;
Smax——中心处的最大沉降值;
i——宽度系数。
最大沉降量采用下式估算:
Vs——盾构引起的地层损失;
i——沉降槽的宽度系数。
计算结果显示,沉降呈正态分布,根据数据作图如5所示。
由云图6可以看出,在盾构周围有较大的应力,为了保证盾构隧道的强度,及时给盾构后的四周支护衬砌。激活衬砌单元,使衬砌的刚度发生作用,目前的地铁施工中,采用拼装管片的方法使圆形管片形成封闭的衬砌环,可以得到后续第二步变形如图7所示。
从云图7可以看出,衬砌的激活对沉降并没有什么影响,衬砌的目的是用了加强隧道壁的刚度的。
按一层楼1.8KN/m2的恒重加载,分别在隧道上方加5层楼和3层楼的荷载,每层楼按100m2折算。
由图9可以看书,加载后,地面沉降有所加强,地铁通常在城市施工,建筑物较为密集,周边建筑物产生的荷载要大的多,为了减小盾构隧道对地面沉降的影响,减小地面建(构)筑物及地下管线受到的破坏,我们对土体进行加固处理,使土层的承载力增大,粘聚力变大,抗剪切强度得到提高,摩擦角增大,实际施工中有采用高压旋喷桩、水泥搅拌桩、高压注浆等方法。
根据图10、11对比未加固土层的云图可以看到,土层经过加固后,变形有了改善,只是在隧道正上方可以看到比较明显的变形,符合peck公式的理论分析结果。
摘要：随着城市经济的发展,地铁等地下工程的施工越来越多,盾构施工使周边建筑物及管线的基础应力得到释放,或由于附加应力的影响引起地基下沉或隆起,这种变形对地面建筑物影响程度与有效间距、施工范围与建筑物间的土质情况等有关。
关键词：地下工程,地基沉降,盾构施工
参考文献
[1]朱伟(译).隧道标准规范(盾构篇)及解说.北京:中国建筑工业出版社2004.196～172.
[2]覃仁辉,隧道工程.重庆:重次大学出版社.2001,1～6.
[3]沈培良,张海波,殷宗泽.上海地区地铁隧道盾构施工地面沉降分析[J].河海大学学报(自然科学版),2003,31.556～559.
[4]张海波,殷宗泽,朱俊高.地铁隧道盾构法施工过程中地层变位的三维有限元分析.岩石力学与工程学报,2005,24期755～760.
车钩钩舌数值仿真分析 篇6
1 不断改进车钩钩舌数值仿真技术的相关措施
1.1 技术人员要熟知有关车钩钩舌数值仿真分析的基本算法原理
为了不断改进车钩钩舌数值仿真技术, 相关的技术人员要熟知有关车钩钩舌数值仿真技术的基本算法原理, 比如有限元静强度分析的基本算法原理。相关的技术人员必须对车钩钩舌的结构强度做好弹塑性分析计算, 要针对具体情况采取不同的分析计算方法。比如对于车钩钩舌结构强度已经给定的应变增量和矢量, 技术人员首先要做到假设这些已经给定的量都是弹性应变增量, 然后应用弹性应变增量对应的计算原理计算出其所对应的应力增量。但技术人员要注意的是在再次检查的时候不能违了屈服准则, 要始终保证结构的应力在屈服面之中和在屈服面上, 从而保证其正常的运作。但技术人员要注意做好增量开始时与结束时应力的变化情况和性质情况, 从而更好地确定其增量中的弹性和弹塑性的成分。除此以外, 技术人员还应该熟知算弹塑性本构关系、塑性应变增量、更新应力和应变以及应力增量等结构受力情况的基本计算原理, 并在此基础上学会车钩钩舌受力分析。
1.2 技术人员要注重对车钩钩舌有限元模型的建模工作
为了做好车钩钩舌有限元模型的建模工作, 技术人员可以借鉴美国EDS公司的I-DEAS Master Series软件系统, 该软件系统具有功能强大的几何造型功能, 相对于以前比较落后的技术而言, 这个建模效率相对是比较高的。技术人员可以借鉴这种比较先进的三维实体造型功能, 同时可以借助其NX I-DEAS创建钩舌实体模型图, 这样比较技术相对比较先进, 安全性比较高, 而且经济成本也比较低。但技术人员在构建这种车钩钩舌模型的时候, 一定不能因为省时省工而忽略车钩钩舌中比较小的圆角和倒角结构的构造, 从而影响整个车钩钩舌的性能。
1.3 技术人员要注重材料属性及评定准则
车钩钩舌采用材料为TB/T2942-1999之C级铸钢, 材料的屈服极限为:σs≥415MPa, 材料的弹性模量E取值为210 000MPa, 泊松比μ取为0.3。并且依据TB/T1335-1996《铁道车辆强度设计及试验鉴定规范》规定, 各载荷工况下, 其应力值 (Von Mises应力) 不应超过材料的屈服极限。
1.4 技术人员要注重计算载荷和约束条件
为了不断改进车钩钩舌数值仿真技术, 相关的技术人员要注重计算载荷和约束条件, 分析并计算好车钩钩舌牵引和压缩载荷数值。以此更好地保证车钩钩舌在受到拉力的时候, 可以保证使其钩舌牵引凸缘承受大部分的拉力, 护销凸缘和钩舌销做好不要不受力。注重在磨合时期尽量保证车钩钩舌的后护销凸缘和钩舌销承受少量的横向力, 但不同部位群殴压力作用方向自然是法向。控制好车钩钩舌承受牵引载荷作用的范围, 这种有两种情况。一种是技术人员将牵引力大小一般控制为1 225k N左右, 其方向为指向钩舌外部平行于钩舌的纵向中心线方向。另外一种是技术人员将车钩钩舌压缩载荷大小控制在1 400k N左右, 其方向为指向钩舌内部平行于钩舌纵向中心线的方向。
1.5技术人员要做好计算结果的分析工作
要不断改进车钩钩舌数值仿真技术, 相关的技术人员要做好计算结果的分析工作, 可以在车钩钩舌承受牵引载荷作用下, 控制好钩舌的Von Mises的应力, 可以将其最大Von Mises应力控制在389MPa左右, 位于钩舌内侧圆弧表面, 小于C级钢材料的屈服极限415MPa, 安全系数为0.94。并将车钩钩舌的剪切应力控制在适当的范围内, 最大剪切应力可以控制在340MPa左右, 位于钩舌内腕面中间部分, 小于C级钢材料的屈服极限415MPa, 安全系数为0.82。在车钩钩舌承受压缩载荷作用下, 技术人员可以将车钩钩舌最大等效应力控制在在钩舌的内腕面两端以及其冲动台曲面上, 并将最大压缩Von Mises应力控制在为400MPa左右, 小于C级钢材料的屈服极限415MPa, 安全系数为0.96。
2结论
不断改进车钩钩舌数值仿真技术, 可以使得车钩在承受牵引载荷和压缩载荷作用下, 车钩钩舌的最大Von Mises应力均小于C级钢材料的屈服极限, 满足TB/T1335-1996《铁道车辆强度设计及试验鉴定规范》及该车设计任务书的要求, 更好地促进我国铁路运输业又好又快的发展。
摘要：为了验证某型车钩钩舌的运用可靠性, 校核车钩钩舌的结构强度, 运用NX I-DEAS三维软件创建车钩钩舌实体模型, 运用Hypermesh软件划分实体网格单元, 运用ANSYS仿真软件的弹塑性分析模块对钩舌进行数值模拟分析, 计算结果表明, 车钩钩舌的Von Mises应力等效应力小于C级钢材料的屈服极限, 满足TB/T1335-1996《铁道车辆强度设计及试验鉴定规范》及该车设计任务书的要求。
关键词：车钩钩舌,结构强度,C级钢,屈服极限
参考文献
基于数值仿真的电价灵敏度分析 篇7
在过去的20多年中, 世界上许多国家开展了以“厂网分开、竞价上网”为主要特征的电力工业市场化改革, 电力企业从垄断型经营模式向竞争性经营模式转变, 电力市场成为国家经济格局中的重要组成部分。其中, 电价作为市场的核心问题, 直接关系到整个市场的发展和稳定。在电力市场环境中, 电价具有较高的波动性和随机性, 其正常波动是电力市场赖以生存和发展的基本条件, 对市场的健康发展具有积极的正面效应;反之, 电价的异常波动将加剧市场成员的运营风险, 并影响电力市场的稳定发展。尤其是2000年美国加州电价的暴涨给美国加州的经济和生活带来了巨大损失, 也引起人们对电价风险问题的普遍关注。然而, 由于电力市场的特殊性, 如何借鉴金融市场风险管理理论进行电价风险分析还有待研究。目前, 已有一些学者对电价问题展开了较为深入的研究, 并取得了丰硕的研究成果, 其研究主要集中在电价预测[1]、波动性分析[2,3,4]、概率分布[5,6,7,8]等方面。
总体而言, 灵敏度分析是常用的风险分析技术, 能够直观地表征自变量与因变量变化之间的内在联系, 这对于电价风险研究而言, 具有重要的理论意义和实用价值。究其原因, 电力市场因素的不确定性具有其客观必然性, 所以, 进行市场因素变化与电价变化之间的定性和定量分析, 以挖掘市场因素不确定性与电价不确定性之间的潜在规律, 实现对电价风险的深层次分析, 具有极其重要的意义。
目前的灵敏度分析方法已在金融资产市场风险分析中得到较为深入的研究和应用[9], 其研究方法大多是直接对目标变量与自变量之间的数学模型进行求导计算灵敏度系数。在电力系统分析领域, 灵敏度分析被广泛应用于系统仿真[10]、保护控制[11,12]、系统稳定性分析[13]等方面, 这些研究大多基于明确的数学模型或系统, 通过灵敏度分析实现对模型变量及参数的动态分析。然而在电价研究方面, 由于价格本身的复杂性和电力市场的特殊性, 对于电价及其影响因素的显示建模存在较大困难, 从而制约了解析法在电价灵敏度分析方面的应用。为解决该问题, 文献[14]提出了电量电价灵敏度矩阵的计算方法, 此方法所求的电价灵敏度是一种静态量度, 对于电价的动态灵敏度问题不适用。
本文从数值仿真的角度, 提出了一种新的电价灵敏度分析模型。该模型首先根据统计学习理论, 建立基于支持向量机的电价仿真模型, 并对电价影响因素进行小扰动模拟, 然后根据智能仿真结果, 计算电价的灵敏度。该模型将支持向量机的黑箱建模技术与灵敏度分析有机结合, 所建模型原理简单、适用性较强。
1模型原理
首先对电价灵敏度模型进行基本定义。设y (t) 为目标交易时段t的电价, x1 (t) , x2 (t) , …, xn (t) 为电价解释变量, 即模型计及的影响因素, 例如历史电价、负荷、网间交换功率、系统可调发电容量等, x (t) =[x1 (t) , x2 (t) , …, xk (t) ]T∈Rk, 且为避免影响因素数据级和单位差异所引起的灵敏度数量级上的差别, x已进行归一化处理。这里定义电价影响因素与电价之间的映射关系为:
简记为y=f (x) 。则电价y的全增量Δy可近似表达为:
式中:βi定义为相对于负荷、历史电价等不同影响因素的电价灵敏度, 且由于x为归一化数据, 所以灵敏度单位与电价单位一致。
进一步将Δy分解为:
于是, 构建电价灵敏度模型的步骤如下:首先, 建立基于支持向量机的电价仿真模型;然后, 利用该模型, 进行电价的数值仿真;最后, 根据仿真样本, 求解灵敏度系数。
2模型实现
针对动态电价灵敏度研究中存在的有限样本和价格本身的复杂性等问题, 本文提出基于支持向量机的黑箱建模技术, 通过对初始状态附近的各市场因素波动及电价响应的数值仿真, 研究影响因素变化与电价变化之间的动态关系, 实现电价灵敏度分析。
2.1建立基于支持向量机的电价回归模型
电力市场中电价与市场因素之间的回归建模属于复杂的高维非线性回归建模问题。在诸多的回归算法中, 支持向量机作为一项新的数据挖掘技术, 借助于最优化方法解决机器学习问题, 有效地克服了传统算法的“维数灾”和“过学习”等问题, 成为解决高维非线性回归问题的新工具, 已在时间序列预测[15,16]方面得到较广泛的应用。因此, 为保证后续灵敏度模型的准确性, 本文提出基于支持向量机的电价仿真模型。设对于任意交易时段t, 给定电价及其影响数据的训练样本集为:
式中:x∈Rk;y∈R。
根据统计学习理论, 模型采用非线性映射Φ将数据映射到高维特征空间, 然后在该空间内进行线性回归。设电价仿真模型为:
采用ε-不敏感函数时, 得到优化问题:
通过引入拉格朗日算子αi, α*i, γi≥0, 得到优化问题的对偶形式[17]:
记K (xi, xj) =ΦT (xi) Φ (xj) 。由于上述优化过程只考虑高维特征空间中的内积运算, 因此可采用核函数, 实现电价非线性回归建模, 建立电价仿真模型如下:
在此模型的基础上, 进一步利用测试集样本进行电价仿真模型的效验。具体地, 基于支持向量机的电价仿真模型实现流程如图1所示, 通过枚举迭代进行寻优, 在参数调整中考察训练误差和测试误差 (见式 (8) ) , 并选取测试精度最高的模型为电价仿真模型。根据该流程, 对各交易时段分别建立模型。
2.2数值仿真和电价灵敏度系数求解
该模块利用前述的电价仿真模型, 实现市场因素变化与电价变化之间的数值仿真。
设市场因素与电价的初始状态为[x (0) , y (0) ], 其中, x (0) =[x1 (0) , x2 (0) , …, xk (0) ]T∈Rk。对于任意的电价影响因素xi (i=1, 2, …, k) , 其第s次迭代新值为xi (s) , 且x
式中:j=1, 2, …, i-1, i+1, …, k。
根据上述样本集绘制电价波动与影响因素xi发生小扰动时的数据曲线, 其切线的斜率即为式 (3) 中电价对市场因素xi的灵敏度系数。以次类推, 建立最终的电价灵敏度模型。
3算例分析
本算例采用2002年10月—12月美国新英格兰现货电能量市场的负荷、网间交换功率和电价实际数据, 构建基于智能仿真的电价灵敏度模型, 并进行有效性验证。其中, 10月和11月的相关数据用于电价仿真模型的训练, 建立12月1日的电价灵敏度模型。
设y为t交易时段的系统边际价格, 电价影响因素x1和x2分别为同时段的负荷和网间交换功率, 则电价灵敏度模型形如Δy=β1Δx1+β2Δx2。首先对负荷和网间交换功率数据进行归一化处理, 根据本文的建模方法, 建立基于支持向量机的电价仿真方程, 其中, 核函数取线性核函数。在此基础上, 对负荷、网间交换功率及电价变化进行智能仿真, 仿真过程中负荷和网间交换功率均采用±3%基值以内之间的均匀波动, 输出为电价波动的仿真结果, 从而建立电价灵敏度模型。所建电价灵敏度模型的相关数据如图2和图3所示。图2为2002年12月1日6:00的电价灵敏度曲线。由图知, 电价与负荷、网间交换功率变化之间存在明显的线性相关性。于是, 求出该交易时段电价对负荷和网间交换功率的灵敏度系数分别为3.972 6美元/ (MW·h) 和-3.312美元/ (MW·h) 。由此发现负荷和网间交换功率对电价变化的不同贡献, 总体而言, 该交易时段电价对负荷和网间交换功率的灵敏度数值大致相同, 但不同因素引起的电价变化却刚好反向。
图3为2002年12月1日24个交易时段的电价灵敏度系数曲线。易见, 各交易时段电价对负荷和网间交换功率的灵敏度系数并不相同。但总体而言, 电价对负荷的灵敏度系数为正数, 说明电价变化与负荷变化方向相同, 即负荷增大时电价有增大的趋势;而在多数交易时段内, 电价对网间交换功率的灵敏度系数却为负数, 且绝对值略小, 说明多数时段内电价变化与网间交换功率变化方向相反, 即网间交换功率增大时电价有减小的趋势。
此外, 算例还将电价对负荷与网间交换功率的灵敏度系数进行了数值比较, 所得的灵敏度系数比率曲线如图4所示。可见, 在绝大多数交易时段, 电价对负荷的灵敏度系数与电价对网间交换功率灵敏度系数的比率都大于1。尤其是在9:00和12:00两个交易时段, 甚至超过10多倍。这说明在这些时段, 负荷对电价的影响占主要地位, 网间交换功率的影响占次要地位。例如12:00的电价灵敏度模型为:
当负荷和网间交换功率都变化3%时, 电价变化0.281 7美元/ (MW·h) , 其中负荷引起的电价变化为0.288 6美元/ (MW·h) , 网间交换功率引起的电价变化为-0.006 9美元/ (MW·h) , 易见2个影响因素对于电价变化贡献的较大差异, 这种波动量分解有助于更好地分析和理解电价波动。另一方面, 当比率远小于1时, 则说明此交易时段网间交换功率占主要地位, 负荷影响占次要地位, 例如在第4个交易时段。
从以上分析可见, 本文灵敏度分析模型原理简单, 易于理解和操作, 且建模过程可解释性强。
4结语
本文基于数值仿真的思想, 提出了一种新的电价灵敏度模型。该模型充分利用支持向量机的黑箱建模技术, 通过对市场因素波动的模拟, 考查不同市场情景下的电价动态响应。再根据激励与响应之间的对应关系, 成功求解出电价灵敏度系数。而且, 本文模型适用于电力市场下的小样本条件, 模型泛化能力较强。但是, 由于本灵敏度分析以回归模型为基础, 电价模型的回归精度会对灵敏度准确性有一定的影响。
仿真数值计算 篇8
严重灾害或战争中大量伤病员主要采用担架后送。担架伤员约占伤员总量的80%以上,且到达救护所需要实施医疗处置的伤员很多,其中需要进行中小处置如清创、更换敷料等处理的占了较大的比例,直接在担架上操作会给急救医护人员带来体位不适,换乘又会使伤员增加二次损伤的几率,而手术或处置室中的手术台或诊疗床装备数量有限。为提高担架伤员救治效率,研究了一种新型担架支架,工作状态如图1所示,装备可直接接驳通用担架,形成现场急救、检伤分类、术前准备、留置观察等直接处置担架伤员的简易医疗救治平台,提高伤员的救治效率。
新型担架支架其主要承力部件由2个X型支架及4支相同的横撑组成,如图2所示,均由铝镁合金管材弯制而成。基于担架的尺寸、操作高度及功能要求,承力部件的长度及折弯角度尺寸已经确定,因此,管材的截面惯性矩即尺寸参数对担架支架的安全可靠性及质量具有决定性的影响。在保证可靠性前提下寻求最小装备质量,确定担架支架管材最适截面尺寸,避免为保证某些设计要求而盲目增加设计裕度,没有现成的经验公式或理论可供直接计算[1]。另外,担架支架的实际结构及载荷的传递情况较为复杂,在数值仿真时载荷及约束无法准确施加。本研究通过变换位移参照系,简化担架支架主要承力部件载荷、约束的施加,并对其进行整体建模。以质量最小化为优化目标,约束装备的最大内应力与最大位移,采用基于易于收敛的序列二次规划(SQP)数值优化设计方法[2,3,4],寻求担架支架管材最适截面尺寸,对优化后的结构进行强度、刚度及稳定性数值计算,与试验结果进行对比,并分析偏差产生的原因。
2 担架支架的约束、载荷分析及理想化
2.1 载荷分析及计算模型的整体建模
担架支架主要承受担架及伤员的垂直向下载荷,要求额定承载能力为1 100 N,加载1 h后担架面下降小于10 mm,安全系数大于3。担架支架主要的承力部件为2个X型支架及支腿部2根横撑,X型支架本身由3个部件组成,将承载载荷传递到地面形成平衡支撑。在这种复杂的载荷及力传递状态下,分别对每一部分进行建模、计算,其结果将难以准确,也无法计算装备整体稳定性。因此,本研究对担架支架主要承力部件2个X型支架及支腿部2根横撑进行整体建模,总载荷分为4个等份施加在上部4个支撑柱上,如图3所示。
2.2 基于变换参照系的约束分析及计算模型简化
以地面为位移参照系,担架支架装配体之间实际约束十分复杂,包括了地面对支腿部2根横撑的部分自由度的约束,以及X型支架关节对2根支撑杆的部分约束,且担架支架关节处的质点对地面也存在位移,导致了需约束处的自由度数量、方向与位置不能确定,无法正确施加约束,数值分析无法进行。采用变换参照系的方法,将位移参照系变为担架支架X型支架关节处,这时支架在关节处的质点位移为0,可方便确定约束状态。在变换参照系的情况下,除去所受额定载荷外,还增加了地面对支架的垂直向上的支撑力,支撑力大小与额定载荷相同。根据以上分析,将关节对支撑杆的约束施加到相应的位置上,建模时可将关节略去。结合2.1一节的分析结果,对模型进行整体参数化建模,根据经验,首先假定优化对象外径尺寸d为32 mm,壁厚h为3 mm。模型及约束、载荷施加如图3所示。
2.3 模型的理想化
由图2主要承力部件可以看出,担架支架的计算模型实际是由受约束的杆系组成的,然而,与最短的支撑管的长度相比,由于管材的直径相对较大,把体系简化为桁架模型将产生错误的计算结果。但是其壁厚参数与最短支撑管尺寸相比相对较小,可将体系理想化为薄壳模型,利用Mechanica中的壳单元进行有限元分析,共形成24个壳对,大大提高计算机计算效率。
3 对担架支架模型的工程分析
3.1 前置分析
前置分析需要在建模、约束及载荷施加后赋予模型材质及进行网格划分。担架支架使用的是某牌号的铝镁合金材料,该材质的泊松比为0.3,杨氏模量为70 GPa。共划分三角形和四边形壳单元网格457个,节点377个。网格划分如图4所示。
3.2 后置处理
在Mechanica中新建静态力学研究项目,采用多通道自适应方式寻找收敛结果,多边形最高逼近阶次为9,计算出Ф32 mm×3 mm模型的应力及结构位移如图5及图6所示。
由数值仿真结果可以看出,最初设计的外径32 mm,壁厚为3 mm时,支架模型的最大内应力为22.6 MPa<[δ]=60MPa,最大内应力发生在长支撑管约束处。在变换坐标系的情况下,支架模型的两端受载荷作用处均有位移,最大位移发生在如图6所示位置,相对于地面位移参照系,模型的最大位移为1.14 mm,此时计算模型的质量为3.42 kg。
4 担架支架的优化设计
4.1 优化设计的数学模型
式中,x1、x2分别为铝镁合金管的外直径以及壁厚;目标函数W(X)为模型的总体质量;σj(X)为结构主要部位的复合应力,共J个;[σ]为材料许用应力;分别为xi的上、下限。
4.2 设计变量的取值范围
设计变量的取值东范围见表1。
mm
4.3 基于SQP的优化设计过程和强度、刚度的试验校核
在Mechanica中建立优化研究,在前述的材质、约束和载荷下,以最小质量为优化目标,同时满足支架最大内应力不大于许用应力60 MPa,确定铝镁合金管截面的最适尺寸。Mechanica计算出铝镁合金管的最佳外径为28.75 mm,最佳壁厚为1.89 mm,依照设计规范进行数据的圆整后,确定管尺寸为Ф30 mm×2 mm,此时支架模型的质量为2.36 kg,最大内应力满足前设条件,其应力及位移的云图如图7和图8所示,其最大内应力为33.5 MPa,依然出现在长支撑管约束处,安全系数为5.37;对地面的最大位移为1.6 mm。随后,对样机进行额定载荷负载试验,结果形变为2 mm,符合设计刚度要求,与计算值1.6 mm接近(见表2)。形成偏差主要原因可能是实际担架支架为多部件装配体,部件之间装配空隙累加,经载荷压实后这些空隙尺寸造成与计算值之间的偏差。
注:变形测试结果2 mm,去掉载荷后残余变形为1 mm;符合设计要求
5 优化后担架支架整体稳定性校核
细长型支撑架在受到轴向压力作用时,虽然强度、刚度符合要求,但可能在远未达到材料屈服强度时出现突然失稳,导致结构破坏而丧失工作能力。因此,对担架支架进行整体稳定性校核是十分必要的。以静力分析为基础,在Mechanica中计算担架支架的三阶失稳模态,并确定机构稳定性安全系数,见表3。模型最低一阶失稳系数为19.5,最低失稳载荷为2 340kg,远大于担架的额定承载,因此在正常工作状态下,担架支架不会出现失稳破坏可能性。图9及图10分别表达了支架的一阶、二阶失稳模态(形变放大100倍),其变形的方向与部位均不相同。
6 结语
本研究针对担架支架工作的典型工况,在变换参照系的情况下,简化了载荷、约束的施加,对担架支架主要承力部件进行了整体建模,以其内置的有限元模块Mechanica为基础,以质量最小化为优化目标,约束装备的内应力与最大位移,通过数值计算方法确定设计变量取值,采用序列二次规划(SQP)数值优化的方法,对担架支架铝镁合金管材截面尺寸进行了优化设计,在满足应力应变约束的条件下,减轻了担架支架的质量。对优化后的结构进行了强度、刚度及稳定性数值仿真校核,均满足设计要求。利用Pro/Mechanica进行优化设计不失为一种很有效的设计方法,可用以完成类似复杂结构的产品设计。
摘要：目的:对新型担架支架的主要结构参数进行优化设计与数值仿真。方法:通过变换位移参照系,简化担架支架主要承力部件载荷和约束的施加,并对其进行整体建模。以质量最小化为优化目标,约束装备的最大内应力与最大位移,采用Pro/Mechanica中序列二次规划(SQP)数值优化设计方法,寻找担架支架管材最适截面尺寸,并对优化后的结构进行强度、刚度及稳定性数值计算。结果:新型担架支架的优化结果能够满足设计要求。结论:利用Pro/Mechanica进行结构优化是一种很有效的设计方法,可为类似机构提供设计指导。
关键词：担架支架,优化设计,数值仿真,Pro/Mechanica
参考文献
[1]宁洁.汽车后送伤员附加装置结构优化设计研究[J].医疗卫生装备,2003,24(9):11-12.
[2]Zhu Zhi-bin.A new SQP method of feasible directions for nonlinear programming[J].Appl Math Comput,2004,148:121-134.
[3]高自友.非线性约束下的SQP可行方法[J].应用数学学报,1995,18(4):579-590.
[4]Proengineer wildfire help,V5.0[CP/CD].Needham MA,USA:PTC,2009.
仿真数值计算 篇9
进入21世纪后, 随着仪器设备和计算机技术的发展, 细观尺度下沥青混合料力学行为的相关研究逐渐成为热点。在构建细观力学数值模型过程中, 有限单元法将沥青混合料视为粗集料和砂浆两相[1], 离散单元法则划分为三类接触:集料内部、砂浆内部、集料与砂浆界面[2]。无论何种数值方法, 首先需获得沥青砂浆的宏观力学参数。通常静态参数如杨氏模量、抗压或弯曲强度等较易直接测得, 而由于试验操作难度、材料本身性质等原因, 砂浆动态参数通常采取对混合料加载过程进行反演获得。
本文借助基于离散单元法的三维颗粒流软件 (Particle Flow Code in Three Dimensional, PFC3D) , 利用室内试验与数值模拟相结合的方法, 对典型级配的沥青砂浆进行动态模量试验, 获取沥青砂浆的动态黏弹性参数, 研究材料高温黏弹性本构和动态力学行为, 同时为沥青混合料细观力学研究提供材料参数。室内试验与数值仿真的相互验证, 还能为后续使用虚拟试验部分代替真实试验提供借鉴。
1 室内试验与数值模型
1.1 沥青砂浆试件
为体现研究成果的普适性, 选取4种矿料级配 (AC-13、AC-20、SMA-10、SMA-13) 作为成型沥青砂浆的依据。文中采用玄武岩和SBS改性沥青作为原材料, 将粒径小于2.36 mm的细集料与沥青结合料的混合物视为砂浆, 根据沥青混合料级配换算出砂浆的各筛孔质量通过率, 即沥青砂浆级配, 如表1所示。
参考Superpave体积设计法2和《公路沥青路面施工技术规范》 (JTG F40—2004) 推荐公式, 假定矿粉与细集料表面均匀裹附着沥青, 颗粒表面积与其所裹附沥青质量为正比例关系, 按照比表面积法[5]获得沥青砂浆的油石比, 如表1所示。
沥青玛蹄脂混合料结构类型属于骨架密实型, 相对密级配混合料而言, 细料比例较少, 但油石比和矿粉含量较高, 因此通过换算得到砂浆的沥青用量偏高。可以认为, SMA系列成型的沥青砂浆的黏弹性能更为显著, 此外, 砂浆较高的沥青用量使得常规试验方法较难直接获取材料的高温黏弹性参数, 需根据时温等效原理予以确定。
文中动态模量试验试件尺寸为高150 mm、直径100 mm的圆柱体, 部分试件采用旋转压实 (SCG) 成型;部分试件在高温下呈现流动状态, 可直接倒入模具, 借助静压法成型, 待冷却后切割为目标尺寸。针对每类级配制备3个平行试件, 共4组12个试件进行测试, 各试件周边粘贴3组位移传感器, 如图1 (a) 所示。根据《公路工程沥青及沥青混合料试验规程》 (JTG E20—2011) 的T0738—2011章节, 通过Superpave简单性能试验机 (SPT) 开展不同温度 (15℃、20℃、25℃及30℃) 与频率 (0.1 Hz、0.5 Hz、1 Hz、5 Hz、10 Hz、25 Hz) 下沥青砂浆的单轴压缩动态模量试验, 如图1 (b) 所示, 对试件施加偏移正弦波轴向压应力, 控制应变模式 (85~110με) 。
1.2 数值模型
在PFC3D内编写程序生成规则排列的单元, 即每单元四周与6个单元相接触, 设定单元半径1 mm, 最终构建如图2 (a) 所示由148 200个离散单元组成的高150 mm、半径100 mm的圆柱体数字试件。文中将砂浆看作均质体, 因此所建立的数值模型由同样大小的单元组成, 并对每一单元赋予相同的微观特性, 相邻单元之间接触本构行为可以看作砂浆单元与砂浆单元之间的接触。
PFC3D软件提供的用户自定义Burger's接触模型如图2 (b) 所示, 由麦克斯韦尔 (Maxwell) 和开尔文 (Kelvin) 两部分模型串联组成, 包含切向和法向作用。Maxwell模型和Kelvin模型被广泛应用在沥青混合料黏弹性分析中, 分别侧重于描述材料的蠕变和应力松弛行为。而沥青结合料/砂浆/混合料的黏弹性能往往呈现出蠕变和松弛的叠合效应, 因此, 综合了Maxwell和Kelvin模型特点的Burger's模型, 结构简单明了、参数数目适宜, 是研究人员最为常用的经典黏弹性模型。
离散单元法通过微观力学响应表征材料宏观力学特征, 室内试验测得的材料参数实际是细观力学行为的宏观反映, 不能直接作为离散元模型的输入参数, 有关Burger's模型中宏观力学属性与微观力学元件之间联系, Liu[6]已在其论文中给出详细推导过程, 文中直接引用, 图2 (b) 所示宏观本构行为与法、切向微观接触行为之间的关系为
式中:E1、E2分别为Maxwell和Kelvin模型内弹簧劲度;η1、η2分别为Maxwell和Kelvin模型内黏壶黏度;L为相邻单元球心距 (文中即为单元直径) , v为沥青砂浆泊松比, μ为摩擦系数;分别为Burger's接触模型法向微观特性;分别为Burger's接触模型切向微观特性。
综上所述, 通过沥青砂浆室内单轴蠕变试验获取宏观Burger's四参数, 借助公式 (1) 和 (2) 即可得到离散元模型的微观特性。
2 结果与分析
2.1 动态模量
PFC3D软件提供了“WALL”命令, 通过控制“WALL”单元的速度实现对模型的加载。在如图2 (a) 所示模型顶底面各生成直径100 mm的圆形“WALL”单元, 固定底面自由度, 不断调整顶面“WALL”单元的移动速度实现正弦波加载。在数值模拟过程中, 实时采集模型顶底面“WALL”单元所受的不平衡反力, 取均值后与作用面积相除即可获得轴向应力;通过实时监控“WALL”单元位移, 与模型原始高度相除即可获得轴向应变。
以AC-13型砂浆为例, 15℃条件下, 0.1 Hz、1 Hz、5 Hz和10 Hz下其应力应变如图3所示。
离散元模拟与室内试验获得的应力应变曲线发展趋势基本一致, 由模拟结果可以看出, 荷载应力按照正弦波形式施加于数值模型;随着频率的增大, 轴向应力幅值逐渐提高。由图3可见轴向应变幅值随着荷载作用时间的延长而呈现出增大趋势 (0.1 Hz低频作用下尚未体现) , 表明与沥青混合料相比, 沥青砂浆具有更为显著的黏弹特性, 在重复荷载作用下的蠕变变形累积效应更为明显。此外, 相比AC型砂浆, SMA型砂浆在相同频率作用下反映出的轴向正弦应力偏小, 此现象也与之前提及SMA型砂浆具有更明显的黏弹特性这一认识相符。
沥青砂浆的动态模量E*和相位角Φ可以通过动态模量室内试验或离散元模拟得到的应力应变曲线并结合式 (3) 求得:
式中:σmax、σmin、εmax和εmin分别为加载应力、应变响应的峰值和谷值;Δt为相邻峰值应力和峰值应变的时间差;T为加载周期, 加载频率的倒数。
基于图3所示应力应变曲线, 按照式 (3) 计算得到离散元模拟的沥青砂浆动态模量和相位角。对15℃、25℃条件下5个加载频率的数据进行比较, 结果如图4所示。
由图4可以看出相同温度与加载频率下, 无论数值模拟还是室内试验, 最大粒径较大的AC-20型砂浆动态模量均略高于AC-13型, 而最大粒径较小的SMA-10型砂浆动态模量却高于SMA-13型;由表1数据可知, SMA-10型砂浆中含有的0.3~0.6 mm、0.6~1.18 mm和1.18~2.36 mm这3档细料分别高于SMA-13型6.5%、5.5%和2.5%, 同样AC-20型砂浆在0.3~2.36 mm区间内的细料含量也略高于AC-13型;此外, AC-20型和SMA-10型砂浆的沥青用量也分别高于AC-13型和SMA-13型。分析认为, 细集料含量的增加仍然提高了沥青砂浆的抵抗变形能力, 而过高的沥青用量增大了砂浆的蠕变变形。
由图4还可以看出, 各温度、不同频率以及不同类型沥青砂浆条件下, 离散元模拟得到的动态模量均高于室内试验结果;分析认为, 虽然砂浆是由粒径2.36mm以下细集料与沥青胶结料组成, 但仍是非均质混合物, 强度存在空间分布的随机性, 而文中数值模型将砂浆作为均质结构考虑, 每一离散单元均赋予相同的微观特性, 是一种理想状态, 因而导致了模拟结果与室内试验的误差, 最大13%, 平均8%。此外, 不同温度条件下, SMA型砂浆模拟与室内结果的匹配程度 (最大误差8%, 最小4.8%, 平均6.2%) 要优于AC型砂浆 (最大误差13%, 最小7.4%, 平均9.5%) , 说明文中引入的Burger's接触本构模型适用于高油石比的间断级配类型沥青砂浆时, 能够更好地发挥本构模型的黏弹性特征。
在未开展材料微观特性试算校核的前提下, 文中离散元模拟结果与室内试验结果匹配度较好, 说明文中建立的数值模型、引入的Burger's接触本构模型以及离散元方法, 能够较好地实现具有黏弹特性的沥青砂浆动态模量数值仿真, 通过不断调整砂浆的微观特性, 可以取得与室内结果更优的匹配度。
2.2主曲线建立
参考AASHTO Design-PP62所述方法, 选用包含于力学经验法的S型函数, 选用AASHTO规范中的二次多项式移位因子。将移位因子方程带入S型函数 (即动态模量主曲线通用方程) , 得到最终形式如下:
式中:α, β, δ, γ, a1, a2均为拟合参数;f为试验温度下加载频率;TR和T分别为参考温度和试验温度。
限于篇幅, 文中选取AC-13和SMA-13型砂浆的室内试验数据, 以20℃为主曲线参考温度, 利用Origin8.5软件根据式 (4) 进行拟合, 相关参数如表2所示。2种不同级配类型的沥青砂浆在各温度条件下的动态模量曲线, 按照移位因子变换后如图5所示。需要说明的是, 在室内试验过程中, 由于沥青用量较高, SMA类型沥青砂浆在温度超出25℃即呈现软化现象, 因此在研究时不考虑30℃。
由图5可以看出, 由于荷载频率范围较小, 沥青砂浆动态模量主曲线没有呈现出与沥青混合料动态模量主曲线形式相似的S形分布, 而是在一定程度上表现出均匀的变化状态。可以认为, 沥青砂浆在结构上确是类似均质特征, 所体现出的宏观力学行为较为稳定, 也为文中在数值模拟时假定砂浆为均质材料提供了支持;而沥青混合料则为典型的非均质结构, 即使处于相同荷载与温度条件下同样呈现出力学行为差异。
3结论
(1) 离散元模拟实现了正弦波加载, 砂浆轴向应变幅值随荷载作用呈现增大趋势, 表明沥青砂浆与混合料相比具有更显著的黏弹特性。
(2) 细集料含量较多的砂浆动态模量较大;离散元模拟的动态模量略高于室内结果, 在不调整材料微观特性前提下, 总体匹配较好;油石比较大的SMA类型砂浆模拟结果与室内试验误差更小。
(3) 沥青砂浆动态模量主曲线一定程度上呈现出均匀变化状态, 而非沥青混合料通常的S形分布, 体现了砂浆的类似均质材料特性;主曲线的建立, 可以突破试验条件限制, 得到沥青砂浆在全频范围内动态模量数值和变化规律。
(4) 引入Burger's接触本构的离散元方法, 能够较好地对重复荷载作用下沥青砂浆的动态模量进行仿真, 并可以进一步推广应用在疲劳试验仿真中。
参考文献
[1]Wang Jiangyang, Qian Zhendong.Indirect tension test of epoxy asphalt mixture using microstructural finite-element model[J].Journal of Southeast University (English Edition) , 2011, 27 (1) :65-69.
[2]王江洋, 钱振东, Wang Lin-bing.环氧沥青混凝土裂纹发展的三维可视化离散元模拟[J].湖南大学学报, 2014, 41 (6) :112-119.
[3]谭凤吉.CRTS II型水泥乳化沥青砂浆灌注质量通病及防治措施[J].铁道建筑, 2012 (4) :139-142.
[4]徐浩, 王平, 徐井芒, 等.反复荷载下水泥乳化沥青砂浆力学性能试验[J].建筑材料学报, 2014, 17 (4) :715-720.
[5]汪海年.沥青混凝土微细观结构及其数值仿真研究[D].西安:长安大学, 2007.
【仿真数值计算】推荐阅读：
有限元数值模拟仿真08-06
耦合数值计算06-09
电场数值计算07-28
数值计算模型08-13
数值计算论文08-29
MATLAB数值计算06-02
CFD数值计算07-22
数值计算方法结课心得07-06
数控加工中的数值计算09-09