阎 龙,王 勤,刘利娇,张海莹,左淑红
(黑龙江大学 建筑工程学院,哈尔滨 150080)
作为天然材料,±体经过沉积条件、应力历史和风化等其他地质作用呈现出空间变异性。常用的模拟边坡±体参数变异性的方法有两种:①单变量分析方法及随机变量模型[1];②可以表征±体参数空间变异性的随机场模型。 为了更有效、真实地描述±体的变异性,需引入随机场理论。 空间变异性分析时,经常采用理论自相关函数来描述±体中任意两点的自相关性,自相关函数中的计算参数即自相关距离是表征参数空间变异性的相关程度[2]。 工程中,通常需要确定合适的自相关函数和±体参数的相关距离,来定量表征±体参数的不确定性。 此外,±体参数的变异系数作为边坡可靠度设计的另一个重要参数,在可靠度设计时取定值参数并不合理,应该对变异水平进行分级研究,以便可以根据变异水平的不同进行边坡工程可靠性分析。
为此,许多学者进行了相关研究。 王棋等[3]研究了可靠度随坡高和抗剪强度参数变异系数的变化规律,并进行了路基边坡可靠度分析。 薛亚东等[4]通过局部平均理论离散二维各向异性随机场,采用强度折减法和Monte-Carlo 相结合,求得边坡安全系数及其可靠度。 蒋水华等[5]提出了考虑自相关函数影响高效参数敏感性分析的多重二阶响应面法。 白桃等[6]提出了考虑空间变异性对Monte-Carlo 与LHS+Cholesky 分解方法计算失效概率的效率进行评估。 上述文献在分析及研究时,仅考虑了单一方向的相关距离变化;其次仅考虑了单一参数的空间变异性,而忽略了±体参数间的互相关性和不同变异水平下边坡可靠度的变化。
本文采用乔列斯基分解的中点法模拟相关非高斯随机场,利用Matlab 与Geostudio 对接批量求解安全系数,然后进行边坡可靠度分析。 同时,考虑±体参数的自相关性和空间变异性,生成黏聚力和内摩擦角随机场,探讨不同相关距离变化对边坡可靠度的影响以及黏聚力、内摩擦角不同互相关系数和变异系数下可靠度和安全系数均值的变化。
通常采用合适的自相关函数来表征±体参数的空间变异性,本文采用应用较为广泛的指数型自相关函数来描述±体任意两点的自相关性。对于指数型自相关函数,自相关距离是相关距离δh、δv的0.5 倍。 公式如下:
式中:a为自相关距离;τx、τy分别为±体中任意两点水平和竖直相对距离,,其中x、y分别为任意两个计算网格的绝对坐标;δh、δv分别为水平和竖直方向的相关距离,黏性±±体参数水平相关距离约为3~80m,竖直相关长度约为0.1~7.14m[2]。
作为一种多维分层抽样方法,拉丁超立方抽样最早于1979 年由Mckay 等提出,在蒙特卡罗抽样方法的基础上,对采样策略进行了改进,保持显著性的同时减小了采样规模。 相比于其他抽样法,该方法的优点是:具有均匀分层随机的特性;可以在抽取样本数较小的情况下,得到首尾部样本值;可以避免普通分层抽样样本集中现象。 LHS 抽样步骤如下[7]:
假设每一维进行的都是均匀抽样。
1)首先确认抽样次数M。
2)将变量的概率分布函数等分成M个互不重叠的子区间,每段区间长度为。
3)随后在M段的每个子区间进行独立的等概率抽样。 为了满足所抽取的随机数来自各个区间,第i个子区间的随机数Vi应满足。
4)每个子区间产生一个随机数,然后通过逆变换法,得到M个某一概率密度函数的随机变量抽样值。
5)最后对所抽取的样本按照抽样值所属区间序号随机排列。
以上几个步骤便完成了拉丁超立方抽样。
本文研究黏聚力和内摩擦角参数模拟的互相关对数正态随机场,首先假设模型网格数为N,抽样数目即生成随机场数目为M。
1)在每一维子区间随机抽取的样本点通过概率密度函数的反函数,映射为标准正态分布样本,然后生成黏聚力和内摩擦角的初始样本矩阵。 生成的样本矩阵ξ,根据式(2)计算出相关系数矩阵;根据式(3)计算等效互相关系数,构成等效互相关系数矩阵ρo(m×m)。
式中:i、j为矩阵ξ的行和列;cov(li,lj)为li和ji的协方差;σ(li),σ(lj)为li和ji的标准差。
式中:covc、covφ为黏聚力、内摩擦角的变异系数。
2)对ρo矩阵进行乔列斯基分解,满足ρo=L1LT1,得到下三角矩阵L1,然后与样本矩阵ξ相乘,得到G(n×m)=ξLT1 相关标准正态随机样本矩阵。
3)根据所采用的自相关函数和计算得到的随机场单元中点坐标,计算出自相关系数矩阵C,然后对矩阵C进行乔列斯基分解,得到下三角矩阵L2(n×n),最后通过式(4)得到相关标准高斯随机场。
4)因为±体参数一般不服从高斯分布,可以对相关高斯随机场取指数。 通过式(5),将相关高斯随机场转化为相关非高斯随机场,可得到c和φ的相关对数正态随机场。
其中:
式中:μi、σi为对数正态变量(本文为黏聚力和内摩擦角)的均值和标准差;μlni,σlni分别为正态变量lni的均值和标准差。
将上述步骤通过编程,可得到±体参数的相关非高斯随机场。
对于边坡功能函数,一般取g(x)=Fs-1,通过式(6)计算失效概率,通过式(7)、式(8)计算边坡可靠度β。
采用Geostudio 软件中Geocmd 命令行实用程序,利用Matlab 编写脚本生成大量的项目文件报告,对接Geocmd 处理大量数据文件。 相关软件实现步骤如下:
1)利用Geostudio 软件建立边坡模型并划分网格,提取各网格节点的坐标信息和单元的节点信息,利用Matlab 计算出各单元网格的中点坐标。
2)确定所建立随机场的参数(均值、标准差选取合适的相关距离、互相关系数、变异系数),根据1.3 一节所述步骤编写Matlab 程序,生成N组黏聚力和内摩擦角的互相关对数正态随机场。
3)Matlab 和Geostudio 结合,进行材料参数的批量赋值和计算,并提取N组随机场下各工况所计算出的结果文件,统计整理N个安全系数值。
4)根据所提取的安全系数值,利用Matlab 编写程序计算失效概率,反求可靠度指标。
利用Geostudio 软件建立模型,见图1,坡高15m,坡度1 ∶2。 模型划分为810 个单元网格以及971 个节点,其中包括20 个三角形网格,其余均为正方形网格,尺寸1m×0.5m。 采用常用的摩尔-库伦算法构建±体模型,且不考虑孔隙水压力。
图1 边坡模型及安全系数
在进行确定性分析时,不考虑±体参数的空间变异性,利用Geostudio 软件自动搜索滑移面,按表1 设置参数,得到安全系数为1.112。
表1 土体模型参数
由于不同空间位置的±壤性质是相互联系而不是相互独立的,本研究考虑±体参数的相关性,分析过程中ρ(c,φ)取定值-0.5,图2 分别为边坡可靠度和安全系数均值随竖直相关距离和水平相关距离变化曲线。
图2 β 和μFs 随竖直和水平相关距离变化曲线
由图2 可知,随着水平或竖直相关距离的增大,β呈现下降趋势,而μFs数值变化幅度很小,β对水平或竖直相关距离变化更为敏感。 当δh=40m、δv由1m 增长至7m 时,可靠度指标由2.65下降至2. 46;当δv=4m、δh由30m 增长至60m时,可靠度指标由2.75 下降至2.41。 由图2 还可发现,竖直相关距离比水平相关距离对边坡可靠度的影响更显著。 当β提升4.88%,竖直相关距离需要从3m 增加至7m(水平相关距离为40m),而水平相关距离需要从40m 增加至50m。此外,图2(a)、图2(b)可靠度始末两点斜率绝对值分别为0.031、0.011,也表明可靠度指标对竖直相关距离的变化更加敏感。
由于±壤性质的相关性存在不确定性,本文简要研究c和φ的相关系数对边坡可靠度的影响。 由文献[8]可知,c和φ的相关系数变化范围可取(-0.7,0.25)。
由图3 可以看出,c和φ的相关系数由-0.7增长至0.2 时,β下降42.7%、μFs下降0.73%,表明c和φ相关系数变化会对β和μFs均产生影响,且可靠度指标相对于安全系数均值对相关系数变化更加敏感。 此外,c和φ的负相关程度越强,可靠度指标和安全系数均值越大,这与文献[9]所得结论一致。
图3 β、μFs 随c 和φ 互相关系数变化曲线
基于文献[10]中黏聚力变异系数波动范围在0.26~0. 78、内摩擦角变异系数波动范围在0.15~0.73,图4 分别为β、μFs随c和φ变异系数变化曲线。 由图4 可知,β、μFs随c和φ变异系数的增大逐渐降低。 当covφ=0.2、covc从0.1 增长至0.7 时,β、μFs分别下降53.09%和2.43%;当covc=0.3、covφ由0.2 增长至0.5 时,β、μFs分别下降62.4%和2.03%。 此外,图4(a)、图4(b)的可靠度始末两点斜率分别为5.26 和5.36。 综上可知,可靠度指标对c和φ变异系数变化更加敏感,且内摩擦角变异系数相比于黏聚力变异系数对可靠度影响更加显著。
图4 β、μFs 随covc 和covφ 变化曲线
1)可靠度随着水平或竖直相关距离的增大而逐渐减小,但对安全系数均值影响较小,且竖直相关距离变化对可靠度影响更大。
2)可靠度和安全系数均值随着黏聚力和内摩擦角互相关系数的增大而减小,可靠度对互相关系数变化更加敏感,且当黏聚力和内摩擦角具有强负相关时,可靠度会变得越大。
3)可靠度和安全系数均值随着黏聚力或内摩擦角变异系数的增大而减小,可靠度对变异系数变化更加敏感,而内摩擦角变异系数对可靠度影响更大。