第6章土的孔隙水压力 6.1概述 研究土坡中的孔隙水压力,首先要弄清在哪些情况下可以通过稳定渗流场的分析来计算 孔隙水压力,哪些情况下不能简单地通过画流网的办法确定孔隙水压力。在5章中已经提到, 土体内的孔隙水压力通常是在下面两种情况下产生的 (1)孔隙水压力是由水的自重形成的渗流场产生的。这一类问题的一个基本特点是土体 的骨架保持不变,因此可以通过稳定或不稳定渗流场的分析计算较好地得到解决。 (2)孔隙水压力是由作用在土体单元上的总应力发生变化导致的。这一种情况仅发生在 压缩性较大、渗透系数较小的土体中。例如,饱和土地基快速开挖或快速填筑,或者均质土 坝库水位骤降的情况。此时土骨架的体积和有效应力都存在着一个从起始状态到新状态过渡 的过程。而粘性土的渗透系数很小,将水挤出,使土的骨架过渡到新的孔隙比,无法在短期 内实现。这样,就可能出现一个随时间消散的附加的孔隙水压力场。这种孔隙水压力,恰是 导致许多工程失事的直接原因。要解决这一类孔隙水压力的问题,则需要引入一些经验或理 论分析方法。此时,一个简单的、偏保守的方法是假定没有任何水排出,在不排水条件下研 究土的孔隙水压力和强度问题。 在上述两种情况下确定孔隙水压的方法构成了在本章详细讨论的主要内容 6.2粘性土的孔隙水压力系数 当作用在饱和粘性土体上的总应力产生增量△G1和△G3时,孔隙水压力的变化可以通过 下式确定( Skempton,1954) △=Bo3+A(△1-△o (61) 这里,假设土体处于三轴状态,即σ2=03。同时,对于饱和土体,B=1。平均有效应 力p'和偏差应力q(或称八面体法向应力和剪应力)为 q=1=0 (63) 如53.1节所述,P的增量Ap'为 y=4y-△M=-(01-40)
第6章 土 的 孔 隙 水 压 力 6. 1 概述 研究土坡中的孔隙水压力 首先要弄清在哪些情况下可以通过稳定渗流场的分析来计算 孔隙水压力 哪些情况下不能简单地通过画流网的办法确定孔隙水压力 在 5 章中已经提到 土体内的孔隙水压力通常是在下面两种情况下产生的 (1) 孔隙水压力是由水的自重形成的渗流场产生的 这一类问题的一个基本特点是土体 的骨架保持不变 因此可以通过稳定或不稳定渗流场的分析计算较好地得到解决 (2) 孔隙水压力是由作用在土体单元上的总应力发生变化导致的 这一种情况仅发生在 压缩性较大 渗透系数较小的土体中 例如 饱和土地基快速开挖或快速填筑 或者均质土 坝库水位骤降的情况 此时土骨架的体积和有效应力都存在着一个从起始状态到新状态过渡 的过程 而粘性土的渗透系数很小 将水挤出 使土的骨架过渡到新的孔隙比 无法在短期 内实现 这样 就可能出现一个随时间消散的附加的孔隙水压力场 这种孔隙水压力 恰是 导致许多工程失事的直接原因 要解决这一类孔隙水压力的问题 则需要引入一些经验或理 论分析方法 此时 一个简单的 偏保守的方法是假定没有任何水排出 在不排水条件下研 究土的孔隙水压力和强度问题 在上述两种情况下确定孔隙水压的方法构成了在本章详细讨论的主要内容 6. 2 粘性土的孔隙水压力系数 当作用在饱和粘性土体上的总应力产生增量∆σ1 和∆σ3 时 孔隙水压力的变化可以通过 下式确定(Skempton, 1954) [ ( (6.1) ∆u = B ∆σ 3 + A ∆σ1 − σ 3 ∆ )] 这里 假设土体处于三轴状态 即σ = 2 σ 3 同时 对于饱和土体 B=1 平均有效应 力 p′和偏差应力 q (或称八面体法向应力和剪应力)为 ( 1 3 2 3 1 p′ = σ ′ + σ ′ ) (6.2) (6.3) q = σ − 1 σ 3 如 5.3.1 节所述 p′的增量∆p′ 为 ( 1 3 ) 3 1 ∆p′ = ∆p − ∆u = ( − A ∆σ ∆− σ ) (6.4)
152土质边坡稳定分析—原理·方法·程序 对于饱和粘性土,如果在应力发生变化时,假定无体积变化,那么,在知道了土的应力 应变关系的条件下,系数A理应通过理论分析得到。例如,假定土体为遵守广义虎克定律的 弹性体,那么土体的体积变化Δ和Ap'之间存在着唯一的对应关系: K.△p 式中:Ks为土骨架的体积模量。由于△=0,可得Ap=0。因此,由式(64)知,纯弹性的土 体 如果土体是遵守某一个相关联流动法则的弹塑性体,那么在Δ=0的前提下,可以推导 出如下表达式(黄文熙,1989) af af Au= Ap+ ap aq a+K (66) 式中:f为屈服函数;A为硬化参数H的函数,具体表达式参见原文 比较式(66)和式(61),可以得到计算孔隙水压力系数A的表达式 K f可f dp aq A A+K|9 ap 土体的孔隙水压力虽然可以通过上述理论途径确定,但是考虑到各种复杂的因素,从工 程实用角度看,系数A和B仍需通过试验来确定(黄文熙,1989年)。 对式61)可作如下变换, △=BAo (68) 其中 B=B[K0+A(1-K0) (69) K0为静止侧压力系数 K0=△G3/△o1 在土石坝中,可以近似地看作△σ1和Δσ3同步增加或减少,△σy△σ1基本保持不变,这样, B可视为常数,其值可通过类似应力途径的室内试验测定 6.3确定孔隙水压力的理论和方法 6.3.1基本方程 土石坝各运用期的孔隙水压力的确定,属于渗流和固结分析的专门问题,这里仅回顾 些基本的概念 在二维问题中,反映流量平衡的微分方程式为
152 土质边坡稳定分析 原理 ⋅ 方法 ⋅ 程序 对于饱和粘性土 如果在应力发生变化时 假定无体积变化 那么 在知道了土的应力 应变关系的条件下 系数 A 理应通过理论分析得到 例如 假定土体为遵守广义虎克定律的 弹性体 那么土体的体积变化∆V 和∆p′ 之间存在着唯一的对应关系 V K V p s∆ ∆ ′ = (6.5) 式中 KS为土骨架的体积模量 由于∆V=0 可得∆p′ =0 因此 由式(6.4)知 纯弹性的土 体 A=1/3 如果土体是遵守某一个相关联流动法则的弹塑性体 那么在∆V=0 的前提下 可以推导 出如下表达式 黄文熙 1989 q p f A K q f p f K u p s s ∆ ∂ ′ ∂ + ∂ ∂ ∂ ′ ∂ ∆ = ∆ + 2 * (6.6) 式中 f 为屈服函数 A* 为硬化参数 H 的函数 具体表达式参见原文 比较式(6.6)和式(6.1) 可以得到计算孔隙水压力系数 A 的表达式 2 * 3 1 ∂ ′ ∂ + ∂ ∂ ∂ ′ ∂ = + p f A K q f p f K A s s (6.7) 土体的孔隙水压力虽然可以通过上述理论途径确定 但是考虑到各种复杂的因素 从工 程实用角度看 系数 A 和 B 仍需通过试验来确定 黄文熙 1989 年) 对式(6.1)可作如下变换 ∆ = B∆σ1 u (6.8) 其中 [ (1 )] B = B K0 + A − K0 (6.9) K0为静止侧压力系数 (6.10) 0 3 1 K = ∆σ ∆/ σ 在土石坝中 可以近似地看作∆σ1和∆σ3同步增加或减少 ∆σ3/∆σ1基本保持不变 这样 B 可视为常数 其值可通过类似应力途径的室内试验测定 6. 3 确定孔隙水压力的理论和方法 6. 3. 1 基本方程 土石坝各运用期的孔隙水压力的确定 属于渗流和固结分析的专门问题 这里仅回顾一 些基本的概念 在二维问题中 反映流量平衡的微分方程式为
拿土的孔隙水压力153 ah. a (Kx)+(K,) I+e at 式中:h为水头;u为孔隙水压力;t为时间;e为孔隙比。 (6)中左边为单位时间流进土体的水量,右边为单位时间该土体的体积变形。 K和K,为x和y方向的渗透系数;y为垂直方向坐标值;γw为水容重;e为孔隙比,e 的变化是由有效应力的增量导致的,而应力增量需要通过求解反映静力平衡的微分方程式获 得。此时的问题,本质上是个固结问题。严格地求解静力和流量平衡,称为比奧理论。如果 对式(61)]引入一些简化条件,可得到以下简化情况 如果骨架的体积压缩模量较大,可以认为不变形,则孔隙水压力主要是由水的自重引起 坝体内的渗流场确定的。这一类情况相应于稳定渗流期,或半透水的砂壳在库水位骤降期。 此时,式(611)右边为零,对于满足达西定律的渗流场,反映流量平衡的微分方程式为 Kx)+(K,x)=0 (6.13) 式(6.13)为稳定渗流或骨架不可压缩土体非稳定渗流的拉普位斯方程。结合相应边界条 件,可用有限元法确定坝体各点的孔隙水压力。这方面,有许多成熟的方法和程序。对于粘 性土,则应考虑式(6.1)的右项,这就是以下第6.32和633节要讨论的内容。 6.3.2太沙基固结理论 1.基本原理 1+e e=m14=m1(A 如果我们假定在荷载变化过程中Ap=0,则有 式中:m为土的压缩系数。这里假定K=K=K,则 Y m, 用微分方程式(615)和初始条件式(6.1)解渗流场,称为太沙基( Terzaghi)理论。这一理论 的核心是假定荷载变化过程中总应力不变 2.太沙基有限元程序FECP2D-1 程序 FECP2D-1是一个可用于计算二维土坝和基础孔隙水压力消散的太沙基固结有限 元程序,它是在二维热传导问题基础上开发的。使用有限元固结程序 FECP2D-1的主要特 点如下
第 6 章 土的孔隙水压力 153 t e y e h K x y h K x x y 1 1 ( ) ( ) ∂ ∂ + = − ∂ ∂ ∂ ∂ + ∂ ∂ ∂ ∂ (6.11) 式中 h 为水头 u 为孔隙水压力 t 为时间 e 为孔隙比 式(6.11)中左边为单位时间流进土体的水量 右边为单位时间该土体的体积变形 y u h w = + γ (6.12) Kx和 Ky为 x 和 y 方向的渗透系数 y 为垂直方向坐标值 w为水容重 e 为孔隙比 e 的变化是由有效应力的增量导致的 而应力增量需要通过求解反映静力平衡的微分方程式获 得 此时的问题 本质上是个固结问题 严格地求解静力和流量平衡 称为比奥理论 如果 对式(6.11)引入一些简化条件 可得到以下简化情况 γ 如果骨架的体积压缩模量较大 可以认为不变形 则孔隙水压力主要是由水的自重引起 坝体内的渗流场确定的 这一类情况相应于稳定渗流期 或半透水的砂壳在库水位骤降期 此时 式(6.11)右边为零 对于满足达西定律的渗流场 反映流量平衡的微分方程式为 ( ) ( ) = 0 ∂ ∂ ∂ ∂ + ∂ ∂ ∂ ∂ y h K x y h K x x y (6.13) 式(6.13)为稳定渗流或骨架不可压缩土体非稳定渗流的拉普位斯方程 结合相应边界条 件 可用有限元法确定坝体各点的孔隙水压力 这方面 有许多成熟的方法和程序 对于粘 性土 则应考虑式(6.11)的右项 这就是以下第 6.3.2 和 6.3.3 节要讨论的内容 6. 3. 2 太沙基固结理论 1. 基本原理 由于 ( ) 1 1 e m p m p u e ∆ = v∆ ′ = v ∆ − ∆ + (6.14) 如果我们假定在荷载变化过程中∆p 0 则有 t u y u x u cv ( ) 2 2 2 2 ∂ ∂ = ∂ ∂ + ∂ ∂ (6.15) 式中 mv为土的压缩系数 这里假定 Kx=Ky=K 则 w v v m K c γ = (6.16) 用微分方程式(6.15)和初始条件式(6.1)解渗流场 称为太沙基(Terzaghi)理论 这一理论 的核心是假定荷载变化过程中总应力不变 2. 太沙基有限元程序 FECP2D−1 程序 FECP2D−1 是一个可用于计算二维土坝和基础孔隙水压力消散的太沙基固结有限 元程序 它是在二维热传导问题基础上开发的 使用有限元固结程序 FECP2D−1 的主要特 点如下
154土质边坡稳定分析—原理·方法·程序 l)施工的任何阶段坝的横断面任何点上的大主总应力σ1变化都能由有限元总应力分析 法求出; 2)对于出现孔隙水压力的坝体部分,固结分析和总应力分析可以使用相同的有限单元 网格; )在分析中,很容易模拟复杂的几何、边界条件和非均匀材料分布; 4)由于计算初始孔隙水压力的总应力取决于有限元应力分析,因而不需要对总应力的 分布作简化假定 程序FECP2D-1由一个主程序和四个子程序组成,主程序和每个子程序逻辑的功能简 单介绍如下。 (1)主程序:调用子程序,读每个施工和停工阶段所需要的数据,并输出计算结果 (2)子程序 GENER:读单元和节点数据,并形成所有的有限元网格,计算方程数和半 带宽。 (3)子程序 ASEMB:形成流动向量,调子程序 ELMAT,获得单元传导矩阵并组装成 总传导矩阵,修正的总流动向量和总传导矩阵使其满足特定的孔隙水压力边界条件 (4)子程序 ELMAT:形成每个单元的传导矩阵,并返回到 ASEMB子程序。本程序使 用具有线性孔隙水压力函数的三角形单元 (5)子程序 SOLVE:采有标准的高斯消去法求解方程,如果几何条件,时间间隔和c 保持不变,利用相同的减化总传导矩阵可求解不同的流动向量的方程,这一过程会节省大量 的计算时间。 下面通过两个例题对程序进行考核。 [例61]与 Koppula and Mogenstern(1972)算例比较。 为了验证程序的可靠性和计算精度,应用如图6.1所示的模型和 Koppula and Mogenstern 的有限差分法进行对比( Eisenstein,1976)。所用参数、初始边界条件见图6.1,所需其它参数 通过前面和下面关系得到: (1+v)(1-2v) E(1-V) E (6.18) 式中:m为一维压缩系数;cr为固结系数;E为关于总应力的杨氏模量:泊松比v=0.49。 不同时段的固结度计算结果和有限差分法对比如图62所示,可见结果精确度较高 [例6,2]与 Gibson(1958)算例比较。 某水坝的平面计算有限元模型如图63所示( Eisenstein,1976),所作假定为:坝高300t°, 分为厚度相等的10层单元。坝体分层连续填筑,速率为10f/月。坝壳料和心墙料具有不同 的弹性模量和相同的容重,y=140b/t3,泊松比v=0.35 为了对比由不同的荷载转移引起的孔压分布,制定了如下4种方案 °本例原文使用的单位为英制,1t=0.305m,lb=4448
154 土质边坡稳定分析 原理 ⋅ 方法 ⋅ 程序 1) 施工的任何阶段坝的横断面任何点上的大主总应力σ1 变化都能由有限元总应力分析 法求出 2) 对于出现孔隙水压力的坝体部分 固结分析和总应力分析可以使用相同的有限单元 网格 3) 在分析中 很容易模拟复杂的几何 边界条件和非均匀材料分布 4) 由于计算初始孔隙水压力的总应力取决于有限元应力分析 因而不需要对总应力的 分布作简化假定 程序 FECP2D−1 由一个主程序和四个子程序组成 主程序和每个子程序逻辑的功能简 单介绍如下 (1) 主程序 调用子程序 读每个施工和停工阶段所需要的数据 并输出计算结果 (2) 子程序 GENER 读单元和节点数据 并形成所有的有限元网格 计算方程数和半 带宽 (3) 子程序 ASEMB 形成流动向量 调子程序 ELMAT 获得单元传导矩阵并组装成 总传导矩阵 修正的总流动向量和总传导矩阵使其满足特定的孔隙水压力边界条件 (4) 子程序 ELMAT 形成每个单元的传导矩阵 并返回到 ASEMB 子程序 本程序使 用具有线性孔隙水压力函数的三角形单元 (5) 子程序 SOLVE 采有标准的高斯消去法求解方程 如果几何条件 时间间隔和 cv 保持不变 利用相同的减化总传导矩阵可求解不同的流动向量的方程 这一过程会节省大量 的计算时间 下面通过两个例题对程序进行考核 [例 6.1] 与 Koppula and Mogenstern (1972)算例比较 为了验证程序的可靠性和计算精度 应用如图 6.1 所示的模型和 Koppula and Mogenstern 的有限差分法进行对比(Eisenstein, 1976) 所用参数 初始边界条件见图 6.1 所需其它参数 通过前面和下面关系得到 (1 ) (1 )(1 2 ) ν ν ν ′ − ′ + ′ − ′ = E mv (6.17) E E′ + ′ − ′ = (1 ) (1 ) ν ν (6.18) 式中 mv为一维压缩系数 cv为固结系数 E 为关于总应力的杨氏模量 泊松比ν′= 0.49 不同时段的固结度计算结果和有限差分法对比如图 6.2 所示 可见结果精确度较高 [例 6.2] 与 Gibson(1958)算例比较 某水坝的平面计算有限元模型如图 6.3 所示(Eisenstein, 1976) 所作假定为 坝高 300ftΟ 分为厚度相等的 10 层单元 坝体分层连续填筑 速率为 10ft/月 坝壳料和心墙料具有不同 的弹性模量和相同的容重 γ =140lb/ft3 泊松比ν′ =0.35 为了对比由不同的荷载转移引起的孔压分布 制定了如下 4 种方案 Ο 本例原文使用的单位为英制 1ft=0.305m, 1lb=4.448N
第6拿土的孔隙水压力155 1)心墙中大主应力等于上覆压力增量即△G1=yM,该方案无需总应力分析; 2)坝壳料和心墙料采用相同的弹性模量,即E壳=E心墙=1000ksf; 3)坝壳的弹性模量是心墙的5倍,即E坝壳5E心墙=5000ks; 4)坝壳的弹性模量是心墙的10倍,即E壳=10E心培=10000ksf。 排水体 粘土心增 排水体:不可压缩 图6.1两边排水的矩形粘土心墙 G1:t时刻平均沉降量 8:总沉降量 有限元解的固结度 差分解的固结度 摘自 Koppula& Morgenstein,1972) 图6.2矩形粘土心墙的国结计算结果 坝壳 图6.3分区坝有限元模型
第 6 章 土的孔隙水压力 155 1) 心墙中大主应力等于上覆压力增量即∆ = σ1 γ∆h 该方案无需总应力分析 2) 坝壳料和心墙料采用相同的弹性模量 即 E 坝壳=E 心墙=1000 ksf 3) 坝壳的弹性模量是心墙的 5 倍 即 E 坝壳=5E 心墙=5000 ksf 4) 坝壳的弹性模量是心墙的 10 倍 即 E 坝壳=10E 心墙=10000 ksf 图 6. 1 两边排水的矩形粘土心墙 图 6. 2 矩形粘土心墙的固结计算结果 图 6. 3 分区坝有限元模型
156土质边坡稳定分析一原理 程序 心墙的固结系数c为50月,孔压系数B为0.5图64给出了不同方案在竣工时中央 心墙的孔压分布图 维孔压消散 m=10fu/mo =30mo ·( gibson解(1958) C.e50ft2/mo 维孔压消散 维孔压消散 E城瓷=E 维孔压消散 E 维孔压消散 孔压比(a/BYH) 图6.4不同方案在竣工时中央心墙的孔压分布图 从图中可以看出,不同的荷载转移对产生坝体内的孔压的影响是明显不同的。计算结果 方案2)和3)差异最大,表明荷载转移在心墙料和坝壳料的杨氏模量较小比率时差异较大。 同时,当E5E心时,心墙中的总应力减小得较慢,因此方案3)和4)的结果差别较小 6.3.3比奥固结理论 1.基本原理 上面提到,太沙基理论的核心是假定荷载变化过程中总应力不变。如果不作这一假定 则式(6.1)右边孔隙比的变化需要引入应力应变关系来确定。这一项是由有效应力的增量导 致的,然而,此值在不知道式(611)左侧的孔隙水压力变化的条件下,是无法确定的。在岩 土工程数值计算中,这类将土的应力应变分析和渗流分析“耦合”的问题称为比奥(Biot)的 理论。通常采用有限元法来求解这类问题,可参阅有关论文( Sandhu and wilson.,1969 Hwang, and morgenstern,1971;陈祖煜,1985)。本节简要介绍其基本理论,在第9章中将介绍采用 有限元求解这类问题的数值方法。有关渗流、固结和应力应变有限元方法的全面论述将在其 它专著中介绍(陈祖煜,2003)。 设某一土体在自重(包括土体中的水重)作用下,在初始时刻=0时处于有效应力场{6} 和孔压4状态。此时土体的位移和应变作为计算参照点,设为零。如果此时受到外荷载的 作用,则土体将产生附加位移,同时还将产生附加的有效应力和孔隙水压力。新建立的孔压 场将改变原有的渗流流速,一部分水从土中挤出。随着时间推移,孔隙水压力和有效应力不 断变化直至最后形成一个与外荷平衡的稳定的{a'}和l。建立一个y坐标向上的坐标系,在 此过程中的任一时刻t,土体应同时满足静力平衡、变形相容和渗流平衡,其微分方程和相 应边界条件描述如下
156 土质边坡稳定分析 原理 ⋅ 方法 ⋅ 程序 心墙的固结系数 cv为 50ft2 /月 孔压系数 B 为 0.5 图 6.4 给出了不同方案在竣工时中央 心墙的孔压分布图 图 6. 4 不同方案在竣工时中央心墙的孔压分布图 从图中可以看出 不同的荷载转移对产生坝体内的孔压的影响是明显不同的 计算结果 方案 2)和 3)差异最大 表明荷载转移在心墙料和坝壳料的杨氏模量较小比率时差异较大 同时 当 E 坝壳>5E 心墙时 心墙中的总应力减小得较慢 因此方案 3)和 4)的结果差别较小 6. 3. 3 比奥固结理论 1. 基本原理 上面提到 太沙基理论的核心是假定荷载变化过程中总应力不变 如果不作这一假定 则式(6.11)右边孔隙比的变化需要引入应力应变关系来确定 这一项是由有效应力的增量导 致的 然而 此值在不知道式(6.11)左侧的孔隙水压力变化的条件下 是无法确定的 在岩 土工程数值计算中 这类将土的应力应变分析和渗流分析 耦合 的问题称为比奥(Biot)的 理论 通常采用有限元法来求解这类问题 可参阅有关论文(Sandhu and Wilson, 1969; Hwang, and Morgenstern, 1971 陈祖煜 1985) 本节简要介绍其基本理论 在第 9 章中将介绍采用 有限元求解这类问题的数值方法 有关渗流 固结和应力应变有限元方法的全面论述将在其 它专著中介绍 陈祖煜 2003 设某一土体在自重 包括土体中的水重 作用下 在初始时刻 t=0 时处于有效应力场{ } 和孔压 u σ 0 ′ 0 状态 此时土体的位移和应变作为计算参照点 设为零 如果此时受到外荷载的 作用 则土体将产生附加位移 同时还将产生附加的有效应力和孔隙水压力 新建立的孔压 场将改变原有的渗流流速 一部分水从土中挤出 随着时间推移 孔隙水压力和有效应力不 断变化直至最后形成一个与外荷平衡的稳定的{σ ′ }和 u 建立一个 y 坐标向上的坐标系 在 此过程中的任一时刻 t 土体应同时满足静力平衡 变形相容和渗流平衡 其微分方程和相 应边界条件描述如下
第6拿土的孔隙水压力157 (1)静力许可的应力场。在弹性力学中已经熟知 ao)=foi (6.19) 其中 U6}=(0 (621) 式中:{G为总应力;{G}为单位体积力;y为土的实际容重,对于饱和土体即为饱和容重 [a]是一个用矩阵表示的偏微分算子符号 ] ax 对饱和土体,总应力和有效应力的关系为 G}={o)+{a (623) 式中:u为孔隙水压力;{a}由下式定义 }2=(110) (624) 上标T表示转置矩阵。 由于达西定律是对水头h而言的,我们把式(623)中的u改为h,根据定义 得知 这里需要注意,在有限元计算中,压应力一概被处理成负的,故具有压力性质的孔压也 被处理成负的 将式(623)、式(626代入式(619),得到 pla}-y回jah-y)={6} (627) (2)变形相容的位移场。位移场Wy=(Wx,W,)代表了物体中任意一点在x和y方向的 位移。当应变很小时,位移和应变满足下式 )=aw (6.28) 式中:{}为应变矢量 (6 位移场的边界条件要求,在所研究物体的一般边界s2上满足下式
第 6 章 土的孔隙水压力 157 (1) 静力许可的应力场 在弹性力学中已经熟知 [ ]{ } { } (6.19) b ∂ σ = f 其中 (6.20) { } ( , , ) x y xy T σ = σ σ τ (6.21) { } ( ) s T bf = 0,−γ 式中 {σ}为总应力 {fb}为单位体积力 γs为土的实际容重 对于饱和土体即为饱和容重 [∂]是一个用矩阵表示的偏微分算子符号 [ ] ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂ = y x x y 0 0 (6.22) 对饱和土体 总应力和有效应力的关系为 { } σ = {σ ′}+{a}u (6.23) 式中 u 为孔隙水压力 {a}由下式定义 { } = (1,1,0) (6.24) T a 上标 T 表示转置矩阵 由于达西定律是对水头 h 而言的 我们把式(6.23)中的 u 改为 h 根据定义 y u h w = − + γ (6.25) 得知 u (h y) (6.26) = −γ w − 这里需要注意 在有限元计算中 压应力一概被处理成负的 故具有压力性质的孔压也 被处理成负的 将式(6.23) 式(6.26)代入式(6.19) 得到 [ ]{ } [ ]{ }( ) { } (6.27) w b ∂ σ ′ −γ ∂ a h − y = f (2) 变形相容的位移场 位移场{ } 代表了物体中任意一点在 x 和 y 方向的 位移 ( , ) x y T W = W W 当应变很小时 位移和应变满足下式 { } [ ] {W} (6.28) T ε = ∂ 式中 {ε}为应变矢量 { } ( , , ) (6.29) x y xy T ε = ε ε γ 位移场的边界条件要求 在所研究物体的一般边界 s2上满足下式
158土质边坡稳定分析一原理·方法,程序 (6.30) 式中:{W}为s2上已知或设定的位移约束条件 (3)本构关系。应力场{)和应变场{}通过本构关系联系起来。对于弹性体, {a'}={△a"+{ob} 其中 △σ}=[C]{e} (632) [C]为根据广义虎克定律建立起来的刚度矩阵,它是对称的。从现在开始,{e}均指由 荷载增量引起的应变增量。 (4)流量连续条件的渗流场。按达西定律,水在土孔隙中的流速矢量}满足下式: }=-K 其中 7=(x,y) (6.34) [K]为渗透系数矩阵,对渗流各向异性的土体 K] k卜 (6.36 K]是对称矩阵,渗流主轴恰为x和y坐标轴 连续条件要求,某一时段内,流出单位土体的水量yy)加上水的压缩量-△等 于土体的压缩量-En=-(Ex+E),可表达为 iV4)+(a T a fE! i au 0 (6.37) o at 式中:Q为不饱和土体的压缩模量,当土体饱和并假定水不可压缩时,1Q为零。 将式(626)式(633代入式(637),得 -(V;(K];h+{a72(s+2mh=0 (638) g at 式(638)是一个含对t微分的方程。求解微分方程通常是在(t0,1)时段内通过积分实现的 因此我们引入卷积的定义。设u和为空间座标x和时间坐标t的函数,则定义卷积( Mikusinski, u*v=u(x, t-r)v(x,r)dr
158 土质边坡稳定分析 原理 ⋅ 方法 ⋅ 程序 {W} = {W} (6.30) 式中 {W}为 s2上已知或设定的位移约束条件 (3) 本构关系 应力场{σ ′}和应变场{ε}通过本构关系联系起来 对于弹性体 { } { } { } (6.31) σ σ σ 0 ′ = ∆ ′ + ′ 其中 {∆σ ′} = [C]{ε} (6.32) [C]为根据广义虎克定律建立起来的刚度矩阵 它是对称的 从现在开始 {ε}均指由 荷载增量引起的应变增量 (4) 流量连续条件的渗流场 按达西定律 水在土孔隙中的流速矢量{V}满足下式 {V} = −[K]{ } ∇ h (6.33) 其中 { } ( , ) (6.34) x y T V = V V ∂ ∂ ∂ ∂ ∇ = x y T { } , (6.35) [K] 为渗透系数矩阵 对渗流各向异性的土体 (6.36) = yx yy xx xy k k k k [K] [K] 是对称矩阵 渗流主轴恰为 x 和 y 坐标轴 连续条件要求 某一时段内 流出单位土体的水量{ } ∇ T {V}∆t 加上水的压缩量 u Q− ∆ 1 等 于土体的压缩量 ( ) v x y −ε = − ε + ε 可表达为 0 1 { } { } { } { } = ∂ ∂ − ∂ ∂ ∇ + t u t Q V a T T ε (6.37) 式中 Q 为不饱和土体的压缩模量 当土体饱和并假定水不可压缩时 1/Q 为零 将式(6.26) 式(6.33)代入式(6.37) 得 { } ([ ]{ } ) { } { } = 0 ∂ ∂ + ∂ ∂ − ∇ ∇ + t h t Q K h a T T γ ω ε (6.38) 式(6.38)是一个含对t 微分的方程 求解微分方程通常是在(t0 ,t) 时段内通过积分实现的 因此我们引入卷积的定义 设u和ν为空间座标x和时间坐标t的函数 则定义卷积(Mikusinski, 1959) (6.39) ∫ ∗ = − t t u u x t x d 0 ν ( , τ )ν ( ,τ ) τ
第6拿土的孔隙水压力159 可得 -Vy;(K]*;V}h)+{a;7{s}+n(h-hn)=0 (5)应力的边界条件。应力场的边界条件要求在所研究物体一段边界s1满足下式 F]=[No)+[Kaju (641) 式中:y=(,元)为作用在s1上已知的单位面积上的力;[a]在式(62)中定义 n 将式626代入式(641),得 行}=[l-Nlan(h-y) (643) (6)渗流场的边界条件。渗流场的边界条件可分为以下三类 第一类边界条件。在所研究的物体边界面s3上,水头满足下式 h=h (64) 式中h为设定的水头值。 第二类边界条件。在所研究的物体边界s4上,流量满足下式 q=[N=-v kiNi (645) 同样,对式(645)在时段(0,1)内积分,得 式中:g=1;q为s4上设定的流量强度 ( n)=(n,, ny) (6 第三类边界条件。所研究的边界为渗流场的自由面 此时,自由面上各点的孔压为大气压,故应满足式(644),其中 (6 式中:y为位置水头的参照高程 同时,在自由面法线方向没有流量,故又应满足式(645),q=0。但是,自由面的位置 通常是未知的,需要在满足这两个条件的情况下予以确定,常需迭代解决。 固结问题的数学提法是,在体力6}和边界上面力{7}的作用下,求满足微分方程式 (627)、式(628)、式(631)、式(640)和边界条件式6.30)式(643)式(644式(646)式(648) 的位移场{四}和水头场h。通常用有限元求解这些偏微分方程组边值问题。在第9章,将简 要讨论具体的求解方法
第 6 章 土的孔隙水压力 159 可得 { } ([ ] { } ) { } { } ( ) 0 − ∇ ∗ ∇ + + h − h0 = Q K h a T T ω γ ε (6.40) (5) 应力的边界条件 应力场的边界条件要求在所研究物体一段边界 s1满足下式 [T ]= [ ] N { } σ ′ + [ ] N { }a u (6.41) 式中 { } ( x y T T = T ,T )为作用在 s1上已知的单位面积上的力 [∂]在式(6.22)中定义 [ ] (6.42) = y x x y n n n n N 0 0 将式(6.26)代入式(6.41) 得 {T } [ ] N { } [ ] N { }a (h y) = σ ′ − γ w − (6.43) (6) 渗流场的边界条件 渗流场的边界条件可分为以下三类 第一类边界条件 在所研究的物体边界面 s3上 水头满足下式 h = h (6.44) 式中h 为设定的水头值 第二类边界条件 在所研究的物体边界 s4上 流量满足下式 q { } V [ ] N { } h[ ] K {n} T T = = − ∇ (6.45) 同样 对式(6.45)在时段(t0 ,t) 内积分 得 g q h K {n} T ∗ = −{∇} ∗[ ] (6.46) 式中 g=1 q 为 s4上设定的流量强度 { } ( , ) (6.47) x y T n = n n 第三类边界条件 所研究的边界为渗流场的自由面 此时 自由面上各点的孔压为大气压 故应满足式(6.44) 其中 0 h = y − y (6.48) 式中 y0为位置水头的参照高程 同时 在自由面法线方向没有流量 故又应满足式(6.45) q = 0 但是 自由面的位置 通常是未知的 需要在满足这两个条件的情况下予以确定 常需迭代解决 固结问题的数学提法是 在体力{fb}和边界上面力{T }的作用下 求满足微分方程式 (6.27) 式(6.28) 式(6.31) 式(6.40)和边界条件式(6.30) 式(6.43) 式(6.44) 式(6.46) 式(6.48) 的位移场{W}和水头场 h 通常用有限元求解这些偏微分方程组边值问题 在第 9 章 将简 要讨论具体的求解方法
160土质边坡稳定分析一原理·方法程序 2.比奧固结分析有限元程序CON2D CON2D是一个用于分析饱和土及非饱和土的应力、位移和孔压的平面应变有限元程序 它求解土体渗流和变形耦合的固结问题,能计算土体在排水、部分排水或完全排水条件下的 变形和孔隙水压力。COND可用于模拟坝体施工期、蓄水期直至稳定渗流情况下的固结过 程以及地基土在回填、建筑物荷载以及储油罐等外力作用下的固结过程。此外,它还能够模 拟开挖过程以及支撑结构的设置和撤除。并可以利用杆单元来模拟加筋土结构以及利用梁单 元来模拟防渗墙等结构 COND程序最初由 Chang和 Duncan于1977年用 Fortran语言编制,其后历经多次 改进:1981年 Duncan, Orazio, Chang, Wong and Nami,1987年 Schaefer and Duncan,1990年 Orazio分别进行了改进工作。 中国水利水电科学研究院对CON2D的源代码进行了大量的分析和测试工作,改正了源 代码中一些产生重要影响的错误(陈祖煜,1985)。这些修正获得了程序编者 Duncan教授本 人的确认,并且在CON2D版本中加入了相应这些修正的源代码,并在注释行中注明系由本 书作者提供。改进后的CON2D能够执行以下功能: 1)进行大规模高效率的有限元计算; 2)提供多种本构模型供实际工程选用; 3)模拟施工、蓄水过程及其后水位涨落变化; 4)对各级荷载按均分比例进行加载 5)按双线法计算浸水湿化变形 6)处理施工及蓄水过程中复杂的边界条件; 7)模拟坝体内部细砂砾石排水层的设置 8)模拟混凝土防渗墙与周围土体的接触特性等 经过多次改进和修正后,CON2D具有如下功能特点 (1)在计算模型及单元类型上,CON2D用于平面应变计算,块体单元类型为四边形等 参单元及退化的三角形单元。单元结点数目可设置为4~8个。接触面单元为两结点单元。 (2)在材料应力应变关系上,材料应力应变关系模型包括:线弹性模型,邓肯一张的E B模型、修正剑桥模型以及刚塑性接触面模型 (3)在非线性有限元方法上,程序采用工程中应用最广泛的中点增量法进行非线性计 算。对每一个荷载増量,程序进行两次计算,利用中点应力状态确定模量矩阵。对每个荷载 增量,还可按比例划分为若干微增量进行加载,使得对施工过程的模拟更为准确。 (4)在存储机制上,程序采用动态内存分配机制,根据计算规模自动确定内存需求量并 进行动态分配。数据交换文件采用二进制格式保证高速读写操作。 3.工程应用实例之 小浪底大坝心墙施工期孔隙水压 小浪底水利枢纽工程位于黄河中游最后一段峡谷的出口处,上距三门峡大坝130km。拦 河大坝为壤土斜心墙堆石坝,设计最大坝髙154m。心墙由粉质壤土组成,防滲墙顶端设有 高塑性士区,坝壳由堆石体组成。河床段坝基处于深厚砂砾石之上。大坝典型剖面详见图 65,建立平面有限元模型如图66
160 土质边坡稳定分析 原理 ⋅ 方法 ⋅ 程序 2. 比奥固结分析有限元程序 CON2D CON2D 是一个用于分析饱和土及非饱和土的应力 位移和孔压的平面应变有限元程序 它求解土体渗流和变形耦合的固结问题 能计算土体在排水 部分排水或完全排水条件下的 变形和孔隙水压力 CON2D 可用于模拟坝体施工期 蓄水期直至稳定渗流情况下的固结过 程以及地基土在回填 建筑物荷载以及储油罐等外力作用下的固结过程 此外 它还能够模 拟开挖过程以及支撑结构的设置和撤除 并可以利用杆单元来模拟加筋土结构以及利用梁单 元来模拟防渗墙等结构 CON2D 程序最初由 Chang 和 Duncan 于 1977 年采用 Fortran 语言编制 其后历经多次 改进 1981 年 Duncan, Orazio, Chang, Wong and Namiq, 1987 年 Schaefer and Duncan, 1990 年 Orazio 分别进行了改进工作 中国水利水电科学研究院对 CON2D 的源代码进行了大量的分析和测试工作 改正了源 代码中一些产生重要影响的错误 陈祖煜 1985 这些修正获得了程序编者 Duncan 教授本 人的确认 并且在 CON2D 版本中加入了相应这些修正的源代码 并在注释行中注明系由本 书作者提供 改进后的 CON2D 能够执行以下功能 1) 进行大规模高效率的有限元计算 2) 提供多种本构模型供实际工程选用 3) 模拟施工 蓄水过程及其后水位涨落变化 4) 对各级荷载按均分比例进行加载 5) 按双线法计算浸水湿化变形 6) 处理施工及蓄水过程中复杂的边界条件 7) 模拟坝体内部细砂砾石排水层的设置 8) 模拟混凝土防渗墙与周围土体的接触特性等 经过多次改进和修正后 CON2D 具有如下功能特点 (1) 在计算模型及单元类型上 CON2D 用于平面应变计算 块体单元类型为四边形等 参单元及退化的三角形单元 单元结点数目可设置为 4 8 个 接触面单元为两结点单元 (2) 在材料应力应变关系上 材料应力应变关系模型包括 线弹性模型 邓肯 张的 E B 模型 修正剑桥模型以及刚塑性接触面模型 (3) 在非线性有限元方法上 程序采用工程中应用最广泛的中点增量法进行非线性计 算 对每一个荷载增量 程序进行两次计算 利用中点应力状态确定模量矩阵 对每个荷载 增量 还可按比例划分为若干微增量进行加载 使得对施工过程的模拟更为准确 (4) 在存储机制上 程序采用动态内存分配机制 根据计算规模自动确定内存需求量并 进行动态分配 数据交换文件采用二进制格式保证高速读写操作 3. 工程应用实例之一——小浪底大坝心墙施工期孔隙水压[25] 小浪底水利枢纽工程位于黄河中游最后一段峡谷的出口处 上距三门峡大坝 130km 拦 河大坝为壤土斜心墙堆石坝 设计最大坝高 154m 心墙由粉质壤土组成 防渗墙顶端设有 高塑性土区 坝壳由堆石体组成 河床段坝基处于深厚砂砾石之上 大坝典型剖面详见图 6.5 建立平面有限元模型如图 6.6