面板堆石坝由于可就地取材、经济、抗震性强等优点，一度成为坝工界青睐的首选坝型[1]，但其结构尺寸相差悬殊，面板最小厚度仅为0.3 m，与坝体整体尺寸差别数百倍至千倍之多。面板作为大坝最重要的结构，应力梯度高，应布置较密的网格[2]，且实际研究[3-4]表明，精细的网格密度将获得更合理的模拟结果。
有限单元法使用简便、通用性强、工程应用广[5]，是面板坝等岩土工程结构安全评价的有效技术手段，但该方法单元形式单一，平面网格严格限制为三角形和四边形单元，网格离散限制多，对复杂边界适应性较差[6]。尤其在尺度跨越大的面板坝工程精细化分析中，需消耗大量时间进行较为繁琐的人机交互前处理，一定程度上降低了自动化程度和分析效率[7-8]。近年来，Wolf等[9]提出了弹性的比例边界多边形有限元(scaled boundary polygon finite element method, SBPFEM)，该方法更加灵活自由，融合了边界元(boundary element method, BEM)和有限元法(finite element method, FEM)的优点，并规避了两者的不足，且已被证明精度和收敛性都高于多边形有限元的数值算法[10]，被学者们应用到多个领域：如电磁问题分析[11]、土与结构相互作用分析[12]、裂纹扩展分析[13-14]、面板动水压力研究[15-16]、三维复杂多面体应用[17-18]、波的传播问题分析[19]等。此外，为探索SBPFEM的非线性应用，Chen等[20]利用半解析的弹性解构造单元形函数，采用两套高斯点方案，实现了比例边界多边形有限元的非线性分析，并将其拓展于多孔介质动力液化分析[21]及三维弹塑性岩土工程应用[22-24]。
采用笔者自主开发的非线性比例边界多边形单元，对典型面板堆石坝进行了系统的静力、动力和永久变形分析计算，并与传统有限单元法对比，以验证提出的方法在面板坝工程分析中的可靠性和合理性。随后，与高效的四分树离散技术无缝耦合，进行面板坝结构快速跨尺度精细化分析应用，以揭示SBPFEM能在面板坝工程精细化分析中脱颖而出的显著优势。
计算程序采用作者完全自主研发的大型岩土工程非线性分析软件平台GEODYNA[25]，非线性比例边界多边形单元被集成到该平台。
比例边界有限元是一种弹性的半解析数值方法，其环向离散、径向解析的特点使其可降低一个计算维度，减少一定的计算量。对任意问题域，可用任意n边形对求解域离散(n > 2)，如图 1所示。通过对每个子域进行求解可获得整个求解域的数值结果。理论详尽推导，可参见文献[9-10, 20]。
弹性框架限制了该方法的应用领域，本着拓展推广SBFEM的原则，采用非线性化的比例边界多边形单元。采用两套高斯点策略：保留原始的边界线高斯点，用于计算相关系数矩阵和半解析的弹性形函数；同时，引入内部高斯积分点，可积分获得考虑了材料非线性的单元刚度矩阵，进而组装得到计算域的总体刚度矩阵，并采用Newton-Raphson迭代算法用于求解非线性平衡方程，从而实现该方法的非线性分析。开发的单元可以求解常规单元及传统常规有限元无法计算的多边形单元(常规凸多边形和四分树单元)，见图 2。详细理论推导和计算可参考文献[9, 20, 26-29]。
采用弹性悬臂梁结构端部受弯分析验证实现程序的正确性和计算精度，结构尺寸及网格如图 3所示，其中l=1.2 m、h=0.4 m、弹性模量E=3×1010、泊松比v=0.2、荷载F=2.4×107 N。
采用SBPFEM和FEM分别计算相同结构，并与理论近似解(采用文献[20]提供的公式)进行对比，分析结果列于表 1。从表 1可以看出，采用的SBPFEM具有更高的精度。
以典型面板堆石坝为例，坝高为155 m，坝顶宽度为17 m，面板坡比1:1.8，主体结构材料分区包括面板、垫层区以及堆石料区。为更合理地模拟面板变形和应力，面板与垫层之间添加无厚度的Goodman接触面单元，面板与趾板之间添加缝单元。常规有限元网格离散单元个数为1 702，节点个数1 816，分析模型见图 4。
静、动力均采用堆石料广义塑性模型[30]模拟筑坝材料变形特性，参数列于表 2。面板与垫层间设置的接触面采用广义塑性接触面本构模型[31]，详细参数见表 3。
面板、趾板及基岩采用线弹性模型，参数见表 4。缝单元参数采用文献[32]建议值，其法向压缩刚度为25 GPa/m，法向拉伸刚度为5 MPa/m，切向刚度为1 MPa/m。静力计算考虑了坝体的施工填筑和蓄水过程，坝体填筑分22个荷载步完成，随后分30个荷载步蓄水至152 m高程。
地震动输入采用规范谱人工波，顺河向峰值加速度为0.2g，竖向峰值加速度为顺河向的2/3。地震波加速度时程见图 5。地震波持续时长为25.00 s，计算时间步间隔Δt取0.02 s。
图 6为满蓄期的坝体位移等值线分布图。从图 6可以看出，FEM计算结果与NSBPFEM结果吻合较好。图 7为满蓄期面板应力变形，FEM与NSBPFEM结果基本一致。
表 5给出了满蓄时坝体和面板应力变形的极值。从表 5可见，两种方法计算结果相差较小，证明NSBPEM在模拟大坝填筑蓄水过程中的精度可以得到保证。
图 8为地震过程中面板顺坡向应力随高程分布规律。FEM计算应力与NSBPFEM计算结果非常接近。
表 6列出了两种方法动力计算中坝顶A点的加速度的极值情况，结果显示两者相差不大，说明了NSBPFEM在动力计算中的可靠性。
图 9为地震后坝体永久变形计算结果。从图 9可以看出，两者变形轮廓基本完全一致，说明NSBPFEM可用于震后永久变形计算分析。
除可用于常规单元(三角形和四边形)分析外，SBPFEM的优势在于其可与传统有限单元不能直接处理的四分树网格离散技术无缝结合，进行高效的跨尺度精细化分析。四分树可快速实现尺度跨越，从1 m级到1 mm级仅需10层递归划分(210 = 1 024)[33]。可自动进行高质量单元离散，且划分网格以正方形为主，仅边界和过渡区为多节点单元。
采用四分树技术，并结合SBPFEM对面板坝进行高效跨尺度精细化分析应用，其中，面板单元尺寸设定为0.5 m，坝体单元大小为8.0 m。如图 10~图 11所示，整个模型剖分时间仅为0.216 s。该技术可高效地实现跨尺度精细化分析，避免了文献[7-8]的不足。值得注意的是，该方法进行网格二次划分时，仅需重新定义最大和最小网格尺寸即可，显著改善了网格离散效率。
静力计算分25个荷载步模拟实际施工过程，分30步缓慢蓄水至152 m高程，静力计算所得的应力变形作为初始应力输入，进行动力时程分析，相应的计算参数列于表 2~表 4。地震动激励输入如图 5所示。
面板为整个堆石坝体系最重要的结构，故其应力-变形情况是工程技术人员最关心的结果之一。图 12给出了满蓄期面板的应力-变形情况，整体挠度较小，最大值出现在1/3坝高附近；应力基本均为压应力，且最大值出现在三分之一坝高附近，此外，还可看出两种方法计算结果较为接近。
地震动作用下，面板动应力分布情况绘制于图 13。从图 13可看出：面板应力有少量差别，其中，SBPFEM计算的最大顺坡向应力最大值为13.2 MPa，FEM计算的最大值则为12.2 MPa。最小顺坡向应力计算中，SBPFEM最大压应力为6.91 MPa，最大拉应力为2.21 MPa；FEM计算的最大压应力为5.92 MPa，最大拉应力为2.96 MPa，详细结果列于表 7。两种方法计算所得规律基本一致，仅数值有少量差别。
联合NSBPFE和四分树网格技术，对面板坝进行跨尺度精细化静动力数值分析，计算结果符合堆石坝正常变形的基本规律，且FEM计算结果与NSBPFE结果很接近，采用跨尺度方案后，对满蓄期坝体位移及面板应力变形影响不大，但对地震作用下面板应力影响较为明显，最大应力差别约1.0 MPa。
总体来看，可以认为NSBPFEM用于跨尺度精细化分析是合理和可行的，可为下一步研究损伤演化和局部破坏提供了技术支撑。
采用自主开发的非线性比例边界多边形单元(NSBPFE)对典型面板堆石坝系统地进行了静力、动力计算分析，结果表明：
1) NSBPFE可简单、便捷地处理三角形、四边形以及传统有限元不能直接处理的五边形、六边形等多边形单元，对堆石坝结构进行系统的静力、动力和永久变形分析，表明该方法结果合理，精度较高，可用于面板坝数值分析。
2) NSBPFE具有更强的灵活性、高效性和通用性优势，主要体现为该方法可与四分树网格技术(多节点单元)自动地无缝结合，以进行面板坝高效的跨尺度精细化分析。
3) 研究工作为进一步分析面板坝防渗面板结构局部损伤演化和渐进破坏提供了有力的技术支撑，且具有较强的通用性，可以拓展到其他土木和水利工程的安全评价。
2019, Vol. 41 Issue (3): 11-18
PDF
2019, Vol. 41 Issue (3): 11-18
PDF