黄土力学参数的相关性对生土窑居可靠度的影响

2013-07-10郭平功童丽萍

河南科技大学学报(自然科学版) 2013年5期

郭平功，童丽萍

(1.郑州大学土木工程学院，河南郑州450001;2.河南城建学院土木工程系，河南平顶山467036)

0 引言

生土窑居是指在原始土中经人工挖凿而形成的穴居居住形式，其支撑体系完全由挖凿成型的纯原状土拱体系作为窑居的自支撑结构，没有栋梁支撑，也没有其他支护，但却能够居而不衰、屹立百年甚至数百年而不坍塌。即使在地震多发区(中国在役窑居大多分布于地震多发地带，45%的窑居区地震烈度在7 度以上)，建造年代百年以上的窑居也很普遍，这充分说明了生土窑居存在之合理、构筑之巧妙，蕴涵着结构和力学奥秘。为全面探讨生土窑居民居的结构性能，充分挖掘生土窑居存在合理性、科学性，对在役生土窑居的可靠度进行分析是非常有效的方法，研究成果将为该类居住结构的安全性评价和寿命估计提供依据。这正是本文的研究意义。

民间在营造生土窑居时，为保证可靠性，十分注意选址。通常选择稳定性较高的离石黄土层，避开黄土的节理面，土层应干燥，排水条件好［1］。其实选址是为了保证黄土材料的强度。黄土的强度由黏聚力c 和内摩擦角φ 提供，两者具有负相关性，通常认为其相关因数的范围为－0.3 ～－0.7［2］。本文讨论c 和φ 的相关性对生土窑居可靠度的影响。

对于含有相关随机变量的结构可靠度问题，早期的一些研究采用正交变换的方法，首先，将相关的随机变量变换为不相关的随机变量，然后用JC 法进行计算。从原理上讲，这种方法是正确的，但计算过于繁琐，特别是需要求矩阵的特征值，不便于应用。近年的一些研究则直接在广义空间(仿射坐标系)内建立求解可靠指标的迭代公式，不需要过多的准备工作，应用简单，是对现有可靠度计算方法的推广［3－4］。

可靠度的计算方法有响应面法、一次二阶矩法、蒙特卡罗法等，生土窑居的功能函数是未知的，适合用响应面法［5］来求解可靠度。响应面法的关键是确定结构的极限状态，本文结合有限元强度折减法［6］，提出“富裕安全因数”(FE)作为生土窑居承载能力极限状态的判别标准。

1 基于有限元强度折减法的FE

1.1 有限元模型

生土窑居的几何参数主要包括窑跨、窑腿宽度、侧墙等，如图1 所示，图1 中A 点为中间窑室的顶点。在营造的过程中，根据窑室使用功能、重要性的不同，窑室分为主窑、角窑、侧窑等，主要的窑跨较大、高度较高，使用空间更大。课题组在河南三门峡陕县进行了实地测量、试验，建立了生土窑居几何参数、黄土力学参数的数据库，其平均值分布如表1、表2 所示。

图1 生土窑居各部位示意图

表1 生土窑居几何参数平均值

表2 黄土力学参数平均值

由文献［7］可知:特征点A 点的竖向位移在按二维和三维分析时很接近，生土窑居可按二维问题研究。为简化问题，窑室几何参数均取平均值，如图1 所示。在建立有限元模型的过程中，需考虑生土窑居的营造过程的独特性，它是在黄土塬或黄土崖上经掏凿而成，并且各窑室的开挖不是同步的。首先，建立窑室未开挖时的有限元模型，只需把有限元模型取得足够大，以忽略边界条件的影响;为减少计算时间，在模型的右边加上对称约束，即只建立对称的一半模型，如图2a 所示，施加重力荷载并平衡初始地应力。然后，开挖中间的窑室，如图2b 所示。最后，开挖两侧的窑室，如图2c 所示。生土窑居中的地坑窑顶面经常有车辆行驶或停留，偏于安全，本文生土窑居的有限元模型均考虑该部分荷载，窑室开挖完成后，在地面施加均布荷载4 kN/m2×1.4 =5.6 kN/m2，其中，1.4 为可变荷载分项因数，标准值取4 kN/m2是参考荷载规范中楼面作为汽车通道及停车库的荷载取值。靠近窑室的部分网格划分的较密，离得越远网格划分越粗。模型底部固定竖向约束，左侧固定水平约束，右侧为对称约束。

本文主要考虑黄土的黏聚力和内摩擦角的相关性的对生土窑居可靠度影响，有限元模型根据计算要求需重复建立多轮。有限元模型的几何尺寸在每轮计算中不需改变，但需改变力学参数。黏聚力c和内摩擦角φ 均为正态随机变量，其标准差分别为5.18 kPa 和2.81°。

图2 有限元分析模型

1.2 基于强度折减法的FE

作为生土窑居承载能力极限状态的判别标准，FE等同于生土窑居的功能函数值，应满足下式:

本文偏于安全地规定“富裕安全因数”与安全因数(以FS表示)关系为:

安全因数为结构或构件的抗破坏强度与设计荷载效应组合的比值，目前，并没有相关的规程或规范对生土窑居的安全因数进行定义。求生土窑居的安全因数是对其稳定性进行分析，稳定性分析是岩土工程中的重要研究领域，有限元强度折减法是研究该问题的有力工具。对于生土窑居的安全因数，实质是黄土的强度储备安全因数。由摩尔－库伦理论，土的抗剪强度由黏聚力和内摩擦角提供，强度折减法即对黏聚力和内摩擦角进行折减，如式(1)所示。

其中，F 为强度折减因数;c、φ 和c'、φ'分别为折减前后的黏聚力和内摩擦角;结构处于承载能力极限状态时，黄土的抗剪强度发挥到最大，此时的折减因数即为安全因数FS［8］。

对于图2 的有限元模型，选取主窑拱圈的中点作为特征点［7］(点A)，当土体的力学参数均取表2 所示的平均值时，特征点竖向位移与F 的关系如表3 所示，位移以向下为正。强度折减因数F 从1.00 逐渐增大，最初每次增大0.10，当竖向位移较大、有可能接近极限状态时，F 每次增大0.01;F 由1.58 增大为1.59 时，向下位移由48 mm 增大为98 mm，发生突变，说明F 为1.58 时窑居处于极限状态，安全因数FS为1.58，FE为0.48。

2 广义随机空间中的响应面法

表3 折减因数与特征点竖向位移关系表

2.1 广义随机空间

随机变量的相关性对结构的可靠度有着明显的影响，通常将相关的随机变量转化为独立的随机变量，如Rosenblatt 变换、正交变换等，但转换过程繁琐，本文直接在广义随机空间内建立可靠指标的迭代公式。广义随机空间和直角坐标系随机变量空间不同之处是，前者坐标轴间的夹角不是直角，而是由随机变量的相关系数确定。求可靠指标与验算点的方法仍采用改进的一次二阶矩法，只需把相关的公式做改变即可。

2.2 响应面方法

生土窑居的功能函数是隐式的，需用响应面函数代表真实的功能函数，然后用改进的一次二阶矩法求其可靠指标与验算点，响应面函数的输入为随机变量c 和φ，输出为FE。通常采用不含交叉项的二次多项式作为响应面函数，该函数具有较高的精度、较少的计算量。分别用X1、X2表示随机变量c、φ 时，响应面函数的形式如式(4)所示。

式中，a、b1、b2、c1、c2为待定系数。响应面方法是一种迭代方法，不是求一次响应面函数就可以得到可靠指标与验算点，需经过多轮有限元计算、建立一系列的响应面函数;每轮有限元计算、建立响应面函数后，使用改进的一次二阶矩法求出该轮的验算点与可靠指标，以得到的验算点为基础，利用式(5)求得下一轮有限元计算的试验中心点。

式中，x(k)、x(k-1)表示第k、k -1 轮的试验中心点;x*(k-1)为第k -1 轮所得出的验算点;第1 轮的试验中心点选为c、φ 的平均值。当前后两轮所得验算点的模满足式(6)时，不再进行有限元计算。

3 生土窑居的可靠度

按照前述方法建立广义空间中的响应面方程，根据式(7)、式(8)和式(9)编制Matlab 程序，求得第k 轮可靠指标β(k)和验算点x*(k)，计算收敛后的可靠指标以β 表示。以c 和φ 的相关系数ρc，φ为－0.3 为例，列出其可靠指标与验算点的求解过程，如表4 所示。文后附第1 轮计算时的Matlab 代码。

由表4 可知:经第1 轮计算，所得验算点为(26.192 6，32.402 2)，对应的富裕安全因数FE为0.02，已接近生土窑居的承载能力极限状态，说明在广义随机空间中，不含交叉项的二次多项式响应面非常准确和高效。用式(5)线性插值得第2 轮计算的中心点(26.109 7，31.558 8)，对应的富裕安全因数FE为0，说明经线性插值后，所得试验点更接近验算点。经两轮计算后，已满足式(6)的要求，不再继续计算。可得c 和φ 的相关系数为－0.3 时，生土窑居的可靠指标为4.375 5，验算点为(26.401 1，31.325 7)。

相关系数ρc，φ分别为－0.7、－0.5、－0.3 和0(即c 和φ 相互独立)时，其可靠指标β、验算点如表5所示。由表5 可知:c 和φ 的负相关程度越大，可靠指标越高。当不考虑两者间的负相关性时，所得可靠指标是保守的。

表4 可靠指标与验算点的求解

当c 和φ 的相关系数不同时，验算点也是不同的，如表5 所示。把表5 中的4 个验算点代入有限元模型中，其对应的FE均为0，说明4 个验算点均正确，生土窑居在不同的γ(φ，c)时均可能处于承载能力极限状态。ρc，φ为－0.7时，验算点与均值点的距离比ρc，φ为其他值时远，验证了负的相关性越大，可靠指标越高的规律。