弹性力学平面问题解法的探讨读后感

弹性力学平面问题解法的探讨读后感

1.弹性力学平面问题解法的探讨读后感 篇一

1 样条虚边界元法

样条虚边界元法是一种间接域外奇点边界元法,是对传统间接奇异边界元法的改进.该法将奇点分布位置由边界移至域外,并采用性态良好的B样条函数去逼近虚荷载函数,提高了间接边界元法的求解精度[1,2,3].

对于弹性力学平面问题,设所研究的弹性域为Ω,其边界为L,域内体力为Fl(l=1,2),如图1所示.将Ω嵌入到无限大域中,在Ω域外作虚边界S,并在其上分布密度未知的虚荷载Xl,根据叠加原理,在真实荷载Fl和虚荷载Xl的共同作用下,无限域内任一点P0处的位移和应力为

式中,Q∈S,Q0∈Ω,和(i,j=1,2;l=1,2)为弹性力学平面问题的Kelvin基本解[4].

将式(1)代入有关的边界条件,得到非奇异虚边界积分方程

式中,P∈L;Hk为边界上的已知函数;为由基本解组成的核函数,它们的具体表达式由不同的边界条件确定;k=1,2表示对于平面问题,边界上存在两个边界条件.

积分方程(2)一般不能获得解析解,它们需要通过数值方法求解.将分布在虚边界S上的未知虚荷载函数展开成三次B样条函数,并采用最小二乘边界子段法消除边界L上的余量,最后建立式(2)的离散化方程为

式中,X为由虚边界S上未知虚荷载样条节点参数组成的列阵;A为X的影响矩阵,其中的元素为核函数与B样条函数乘积的二维积分,可以用高斯数值积分求出;B和C为已知列阵,分别取决于Ω内的体力和L上的边界条件.

求解式(3)得出X后,将其代入式(1)的离散形式,即可求出Ω内任一点的位移和应力.

由于样条虚边界元法的基本未知量仅出现在虚边界上,因此,与需要全域离散的有限元法相比,样条虚边界元法具有建模简单,未知量少,计算效率高等优点;而与传统边界元法相比,它避开了边界奇点所带来的奇异积分、边界层效应等复杂问题,提高了边界及其附近区域上解的精度.因此,在结构可靠度分析中,采用样条虚边界元法作为蒙特卡罗法的样本试验方法具有一定的计算优势.

2 Taylor-Neumann展开重要抽样蒙特卡罗样条虚边界元法

蒙特卡罗法通过随机模拟和统计试验来进行结构可靠度分析,回避了可靠度分析中的数学困难,不受随机变量分布形式与功能函数形式复杂性的影响,应用范围相当广泛.该法的缺点是计算量大,当缺乏功能函数显式表达式而需借助数值方法(如有限元法、边界元法等)求解时,计算量大的缺点更为突出.为此,本文在蒙特卡罗边界元法中引入Taylor展开和Neumann展开技术,以降低单次样本计算时间,同时引入重要抽样技术以减少抽取样本数.这些改进方法的引入大大降低了蒙特卡罗法所耗费的计算时间,并能保持良好的计算精度.

2.1 采用Taylor展开技术计算影响矩阵

在样条虚边界元法的求解过程中,影响矩阵A的形成较为复杂,需要进行大量数值积分,耗费时间相对较多.因此,本文采用Taylor展开技术来形成各组样本下的影响矩阵A,以缩短单次样本试验中A的计算时间.

在第n次样本计算中,记式(3)中的影响矩阵A为An,其二阶Taylor展开式如下

式中,N为抽取样本总数;M为随机变量总个数;ri为第i个随机变量(i=1,2,…,M),rin为第n组样本中第i个随机变量的样本值(n=1,2,…,N),ri0为第i个随机变量的均值,r0=[r1,0 r2,0…rM0]T;A0,和分别为当各随机变量取为均值时的影响矩阵、一阶导数矩阵和二阶导数矩阵.

由上可见,只要计算一次A0,和,即可由式(4)通过简单的矩阵数乘和矩阵相加运算计算各组样本下的影响矩阵An(n=1,2,…,N),这比通过大量数值积分直接计算An快得多.

2.2 采用Neumann展开技术计算影响矩阵的逆矩阵

虽然样条虚边界元法中影响矩阵A的规模较小,但很多时候A为满秩矩阵,所以其求逆过程在样条虚边界元法计算中也比较耗时.为此,本文参照Neumann展开蒙特卡罗有限元法的处理方式[5,6],引入Neumann展开技术解决矩阵求逆的效率问题.具体做法如下:

首先对均值影响矩阵A0求逆,得到其逆矩阵,然后对每次样本试验产生的影响矩阵An,令

则含前m项的展开式如下

由式(3),并利用上式,可得第n组样本下样条节点参数列阵Xn的前m项展开式解答为

式中

展开式的收敛标准可取为[6]

式中‖·‖2表示向量范数,定义为

误差限ε一般取0.01～0.05.

从上述公式可以看出,引入Neumann级数展开式后,在反复求解样条虚边界元法基本方程式(3)时只需对均值影响矩阵A0求逆一次.可以预计,Neumann级数展开法将可有效地提高式(3)的求解效率.

2.3 采用重要抽样技术减少样本数

对于前述M个随机变量ri(i=1,2,…,M),设其联合概率密度函数为f(r1,r2,…,rM),对应结构某一极限状态的功能函数为Z=G(r1,r2,…,rM),则结构的失效概率可表示为

式中D为随机变量的定义区域,I[G(r1,r2,…,rM)]为示性函数

假若存在一个重要抽样密度函数h(r1,r2,…,rM),则式(11)可以写成重要抽样形式为

若以h(r1,r2,…,rM)对(r1,r2,…,rM)进行抽样,则由上式可知失效概率Pf的无偏估计值为

式中,N为抽取样本总数,(r1n,r2n,…,rMn)为第n组样本值(n=1,2,…,N).

为达到降低抽取样本数的目的,选择的重要抽样密度函数h(r1,r2,…,rM)应增大抽取样本点落入失效域的机会,但过多的样本点落入失效域会使计算结果误差增大,为了在提高效率的同时保持计算精度,通常会将重要抽样的重心选择在设计验算点处[7,8].本文先通过响应面法[9]计算得到设计验算点,然后以此作为抽样重心进行重要抽样.具体来说,进行重要抽样的随机变量取为相互独立随机变量,它们的分布类型取为对应的原随机变量的分布类型,方差取为对应的原随机变量的方差,而均值则取为设计验算点所对应的坐标值(即均值点取为设计验算点),据此即可确定重要抽样密度函数h(r1,r2,…,rM)的表达式.

2.4 计算步骤

综上所述,Taylor-Neumann展开重要抽样蒙特卡罗样条虚边界元法的计算步骤如下:

(1)采用响应面法求出设计验算点的坐标值.

(2)取设计验算点为重要抽样的均值点,将其坐标值代入样条虚边界元法列式,形成均值影响矩阵A0以及均值值导数矩陈和,并计算逆矩阵.

(3)以设计验算点作为抽样重心,形成重要抽样密度函数,并进行N次重要抽样,得到N组样本.

(4)由式(4)形成第n组样本下的影响矩阵An(n=1,2,…,N),并利用式(7)计算该组样本下的Xn,将其代回式(1)的离散形式后即可求得该组样本下结构的响应值.

(5)将各组样本下的结构响应值代入功能函数,然后利用式(14)计算结构的失效概率.

3 算例分析

3.1 含中心圆孔方形板受单向均布拉力作用

含中心圆孔的方形板受到单向均布拉力的作用,板宽2H=50cm,孔半径r=5cm,如图2所示.由于问题的对称性,取其1/4部分进行研究,如图3所示.考虑板的弹性模量E泊松比μ,厚度t以及拉力q为相互独立的随机变量,其统计特性如表1所示.该板的强度极限值[σ]=127N/Cm2,当圆孔上侧边缘点D处(见图2和图3)的σx超过[σ]时方形板失效.采用蒙特卡罗样条虚边界元法计算该方形板的失效概率.在样条虚边界元法计算中,取98个边界子段和28个虚边界元,数值积分时采用单向4个高斯点,虚边界与真实边界之间的距离为1 cm.

考虑5种计算工况,各工况考虑的因素如表2所示,计算结果和计算时间列于表3中.

由表3可以看出,在蒙特卡罗样条虚边界元法中分别引入Taylor展开、Neumann展开和重要抽样技术,其计算结果均能保持良好的精度,且计算效率均有不同程度的提高.其中,Taylor展开技术对提高计算效率的作用最为显著,引入后计算时间缩短至工况1计算时间的2.0%;重要抽样技术也发挥了相当好的效果,引入后计算时间缩短为工况1计算时间的19.5%;Neumann展开技术的作用稍弱,仅将计算时间缩短了9.2%.综合采用这3种改进技术,计算效果最好,可以在保持良好计算精度的同时,将计算时间缩短至工况1计算时间的0.7%,计算效率大大提高.

下面考察随机变量变异系数对Taylor展开计算精度的影响.取μ=0.2,t=1.0cm,q=30N/cm,仅考虑E的变异,其均值取为3.15×106N/cm2.当E的变异系数δ取不同数值时,工况1和工况2计算得到的点D处竖向位移v的均值和标准差列于表4中.

由表4可见,工况2的计算误差随着变异系数的增大而增大,这是由于工况2采用了Taylor展开截断式所致.当变异系数增大时,式中被舍去的高阶量影响增大,导致工况2的计算误差增大.但是,由于本文方法采用了二阶Taylor展开式,因此在变异系数达到0.25的情况下,仍然保持了良好的计算精度.

3.2 固支深梁受均布力作用

固支深梁受到均布力的作用,如图4所示.考虑板的弹性模量E、泊松比μ、厚度t以及均布线荷载q为相互独立的随机变量,其统计特性如表5所示.当深梁上边界中点D处(见图4)的竖向位移v超过极限值[v]=0.0228时深梁失效.采用蒙特卡罗样条虚边界元法计算该深梁的失效概率.在样条虚边界元法计算中,取20个边界子段和8个虚边界元,数值积分时采用单向3个高斯点,虚边界与真实边界之间的距离为2.

考虑2种计算工况,各工况考虑的因素如表6所示,计算结果和计算时间列于表7中.

由表7可以看出,本算例的失效概率较低,这种情况下本文方法计算效率提高幅度更大,工况2的计算时间仅为工况1的0.35%.此外,本文方法在用响应面法求出设计验算点的同时,实际上也可以结合一次可靠度方法得到可靠指标β的值.但当功能函数非线性程度较高时,该值将会出现较大误差.对于本算例,由响应面法直接求得的β=3.389,与相对精确解(工况1的结果)反算得到的β=3.033相比,误差达到11.74%.

4 结论

当采用蒙特卡罗样条虚边界元法进行弹性力学平面问题可靠度分析时,可以利用Taylor展开式,由均值矩阵生成不同样本下样条虚边界元法的影响矩阵,避免了大量的数值积分运算;可以引入Neumann展开技术,解决大量样本下影响矩阵的求逆问题;还可以采用重要抽样技术,进一步减少抽取样本的数量.这3种方法的综合运用,大大提高了蒙特卡罗样条虚边界元法的计算效率.此外,由于在Taylor展开式中保留了二阶项,因此,本文方法在结构参数中、小变异情况下均可保持良好的计算精度.

参考文献

[1] Su C,Han DJ.Multidomain SFBEM and its application in elastic plane problems.Journal of Engineering Mechanics, ASCE,2000,126(10):1057-1063

[2]苏成,韩大建.高层建筑筏板基础分析的样条虚边界元法.土木工程学报,2001,34(1):61-66(Su Cheng,Han Dajian. Analysis of raft foundation by spline fictitious boundary element method.China Civil Engineering Journal,2001, 34(1):61-66(in Chinese))

[3]苏成,郑淳.基于Erdogan基本解边界元法计算应力强度因子.力学学报,2007,39(1):93-99(Su Cheng,Zheng Chun. Calculation of stress intensity factors by boundary element method based on Erdogan fundamental solutions.Chinese Journal of Theoretical and Applied Mechanics,2007, 39(1):93-99(in Chinese))

[4] Banerjee PK,Butterfield R.Boundary Element Method in Engineering Science.New York:McGraw-Hill,1981

[5] Shinozuka M,Deodatis G.Response variability of stochastic finite element systems.Journal of Engineering Mechanics, ASCE,1988,114(3):499-519

[6] Yamazaki F,Shinozuka M,Dasgupta G.Neumann expansion for stochastic finite element systems.Journal of Engineering Mechanics,ASCE,1988,114(8):1335-1354

[7]金伟良.结构可靠度数值模拟的新方法.建筑结构学报,1996, 14(3):63-72(Jin Weiliang.A new approach to numerical simulation of structural reliability analysis.Journal of Building Structures,1996,14(3):63-72(in Chinese))

[8]贡金鑫,仲伟秋,赵国藩.工程结构可靠性基本理论的发展与应用(1).建筑结构学报,2002,23(4):2-9(Gong Jinxin,Zhong Weiqiu,Zhao Guofan.Developments and applications of reliability theories for engineering structures(1).Journal of Building Structures,2002,23(4):2-9(in Chinese))