FLAC3D动力分析中的人工透射边界和地震波施加方法

FLAC3D动力分析中的人工透射边界和地震波施加方法从动力学的角度上看，动力响应是确定惯性（质量效应）和阻尼起着重要作用时质点或质点系动力学特性和响应的技术，它包括自振、冲击、谐振动、随机振动等分支。

动力学最早应用于结构抗震设计，自上世纪50年代逐步借鉴到岩土抗震设计中。

动力发展历程可总结为静力理论，反应谱理论和时程分析理论三个阶段。

我们知道，地震的三要素为振幅、频谱和持时。

静力理论只考虑了地震引起的最大振幅，属于拟静力法；反应谱理论考虑了振幅和频谱，但在设计中仍然把地震惯性力视为静力，只能算准动力法；时程分析理论考虑了振幅、频谱和持时，是严格意义上的动力分析法。

通常时程动力分析选用的地震波来自：（1）根据设计反应谱人工合成的场地波；（2）场地附近地震台记录的实测地震波。

由于实测地震波中掺杂了许多噪声和干扰信号，因此在使用前必须滤波去噪、频谱分析、积分变换和基线修正。

滤波去噪是为了消除噪声和高频波，频谱分析是为了检测地震波持时内所含的频率分量和振幅，积分变换可以转换地震加速度波为速度波或位移波，基线修正则是为了消除非平稳地震波中的弹性位移零线漂移、基线偏移等现象，大崎顺彦在其著作《地震动的谱分析入门》中做了详细而生动的说明，并附出了地震波处理的Fortran源程序。

鉴于FLAC3D软件是岩土领域广泛应用的时程动力分析软件，这里以著名的埃尔森特罗波（El Centro）为输入激励，研究基于FLAC3D软件的地震波处理和计算方法。

网站“http://www. /data.htm”提供了31秒的El Centro加速度波数据。

有兴趣者可按《地震动的谱分析入门》的方法选取了前8秒的地震加速度波（共401个记录），然后补零配成了512个记录的加速度波以采用快速傅里叶变换法，首先采用FLAC3D Fish函数库的filter函数进行滤波去噪，然后采用fft函数进行快速傅里叶变换，得到傅里叶加速度谱和功率谱，接着采用integrate函数积分两次求得速度波和位移波，并计算地震位移零线漂移值。

具体可细查flac中的帮助和fish命令流。

由于频谱分析，看似混乱无章的地震波也可以分解成不同频率不同振幅的简谐波的组合，因此，只要我们把握了最简单的简谐波动力响应，对地震波时程响应分析也一目了然了。

动力计算是较复杂的力学分析过程，影响因素包括：输入激励（频谱、振幅、持时）、岩土参数、本构模型、透射边界、模型网格、求解方法等。

直接进行一个复杂边坡的动力计算，由于影响计算结果的因素较多，不利于工程师把握本质的规律，以致无法判定结果的适宜性。

对于大型边坡工程来说，抗震分析需要引起足够的重视，每一步动力计算和响应模拟都必须经得起考证。

因此，本文首先从概念模型出发（简单边坡模型），定量分析边坡的动力响应规律，采用理论解进行验证，考证动力计算过程的正确性。

在此基础上，再计入更复杂的地质条件、边界条件和实际边坡模型，最终达到合理可靠地应用于实际边坡工程中的愿景。

此主题相关图片如下：边坡概念模型.jpg可以从四种典型的概念边坡模型来理解边坡动力响应，①代表退化为一竖直柱体边坡；②代表发育有一个软弱夹层的竖直柱体边坡；③代表发育有一软弱夹层的边坡；④代表发育两组正交优势节理的边坡。

①、②模型中坡顶仅一个自由面，是最简单的概念模型；①、②、③模型都是考虑少数结构面的连续介质模型，是目前数值方法中最常用的形式，本文拟采用Flac3D软件模拟；④模型考虑了较多的结构面，常规有限元法等无法模拟，则根据Edelbro（2003）的建议，采用离散元法（3DEC)更为合理。

关于地震波的选择，①、②模型采用简谐波，以方便对比理论解；③、④模型采用实测地震波（如El Centro波）。

此主题相关图片如下：数值方法的选取.jpg数值方法的选取，修改自Edelbro(2003)简谐波或地震波的施加，一般来说可分为加速度波、速度波和应力波。

这三种输入方法是等价的，采取哪种方法视岩土介质类型和动力边界条件而定。

边坡基底岩土类型一般可分为刚性地基和柔性地基。

如果基底为基岩或变形模量较大的岩土体，可以在底部直接施加加速度波或速度波。

若基底为土体，尤其是软土，FLAC3d帮助手册中建议采用应力波更为合适。

下面以一个简单的正弦加速度波来推导应力波的换算方法，然后笔算出应力波在①模型中的传播时间及产生的质点响应速度和剪应力大小，然后用FLAC3D数值解进行验证。

设正弦加速度波的周期为T、地震持时为t，加速度波为a(t)=(π/T)*Sin(2πt/T) ，则积分后可得速度波v(t)=0.5(1-cos(2πt/T),将速度波换算为压缩应力波σp(t)和剪切应力波σs(t)时，Lysmer和Kuhlemeyer （1969）的公式，则有σs(t)=-ρ*Cs*Vs=0.5(cos(2πt/T-1)*sqrt[Eρ/2(1+μ)]。

E为岩土动弹模，ν为岩土动泊松比，ρ为岩土密度，Cs为岩土介质的S波（剪切波）波速。

假设①边坡的高度H为50m，边坡动弹模E为180/7MPa，动泊松比υ为2/7，密度ρ为1000kg/m3，则剪切波在边坡中的传播速度为100m/s，剪切波从边坡底部到达边坡顶部的时间t为0.5s。

无反射波叠加情形下，边坡的最大水平速度为1.0m/s。

边坡遭遇的最大剪应力为0.1MPa。

采用flac3d计算①模型表明：0.25s时，完成一个周期剪切波的施加，随后，剪切波竖直向上传播，0.5s 时到达顶部，并在坡顶产生自由面放大效应和反射，即在0.5s~0.75s间蓝线的最大值约2.0m/s。

剪切波返回模型底部后，由于没有施加人工透射边界，1.25s后剪切波在底部边界反射回顶部，之后反复振荡。

分析边坡剪应力分布规律可知，最大剪应力为0.1MPa，其位置随着剪切波移动而移动。

将加速度波换算为应力波，边坡底部设置了粘滞吸收边界，分析可知，在应力波反射回底部之前，与加速度波规律一致，不存在剪切波反射，1.25s后应力波被粘滞吸收边界吸收，不产生反射，剪应力衰减为0。

综上所述，采用Fish 命令流模拟的剪切波传播和反射规律与理论解一致，边坡最大剪应力分布规律与理论解吻合，剪切波在边坡自由面上产生放大效应，在约束边界产生地震波反射，采用粘滞吸收边界则消除了地震波反射现象。

此主题相关图片如下：无人工透射边界.bmp无人工透射边界，采用速度波施加此主题相关图片如下：黏滞吸收边界.bmp设置黏滞吸收边界，采用应力波施加newconf dyngrid 1 1 50gen zone copy 15 0 0;;def setupomega = 2.0 * pi * freqpulse = 1.0 / freqendset freq=4.0setup;def waveif dytime > pulsewave = 0.0elsewave = 0.5 * (1.0 - cos(omega * dytime)) endifend;def dyn_timearray vec(3)vec(1) = 18.0vec(2) = 0.0vec(3) = 5.0oo = set_fontsize(1.0)oo = draw_string(vec,'Time')vec(3) = 0.0oo = draw_string(vec,string(dytime))vec(3) = 45.0oo = set_fontsize(1.0)oo = draw_string(vec,'stress')vec(1) = 3.0oo = draw_string(vec,'vel.')endset dyn dt=2e-4;;;;;材料本构mod elasprop shear 1e7 bulk 2e7 dens 1000;range name bottom z -0.1 0.1fix z range z .5 55;;;;;速度时程输入apply nquiet dquiet squiet range bottomapply xvel 1.0 hist wave yvel 0.0 zvel 0.0 range bottom;;----------------------------------------------------------------------;;;施加粘滞边界(必须在自由场边界之后，否则失效)，必须于应力或力时程合用。

;;apply dquiet squiet range bottom;;;;;等价于加速度时程输入;apply sxz -2.0e5 hist wave syz 0.0 szz 0.0 range bottom;;----------------------------------------------------------------------;;-------------------- 设置出图格式his gp xvel 0 0 0his gp xvel 0 0 25his gp xvel 0 0 50his dytimeplo crea haitangplo add surf redplo add vel max 0.5 range x -1 2plo add bcon sxz int 1e4 max 1e5 range x 5 18plo add fish dyn_time blackplo set per offplo set bac whiteplo add hist 1 2 3 vs 4plo showTitle ' 正弦波在无阻尼弹性岩基中的传播'set movie avi size 2000 1000 frameperiod 100 file 01_ela_vel_quiet size 1024 768 step 50 movie start;;------------------ 设置动力计算时间solve age 1.5movi finish。