二维弧形一致黏弹性边界单元适用性问题研究

叶 丹，周建庭，尹尚之

(1. 重庆第二师范学院 旅游与服务管理学院，重庆400065； 2. 重庆交通大学 土木建筑学院，重庆 400074； 3. 中建五局第三建设有限公司中南分公司，湖南 长沙 410000)

0 引 言

在大型结构抗震分析中，半无限域的截断面处理方式是决定土体-结构相互作用计算精度的关键，且设置人工边界条件常用于半无限域空间模型。J.LYSMER等[1]通过在截断面处布置阻尼，完成了黏性人工边界的施加，其概念清晰、操作简便，但存在低频失稳和精度不高等缺陷。为克服以上问题，刘晶波等[2-3]基于柱面波动方程，通过在边界处施加弹簧-阻尼器，建立了具有良好弹性恢复性能的黏弹性人工边界(VSBE)，克服了低频失稳问题，获得了良好的稳定性，目前已在有限元软件中得到了广泛应用。但是，在解决复杂二维半无限域动力问题时，由于大部分有限元软件并没有直接施加黏弹性边界的版块，研究者需要逐个进行弹簧-阻尼器的计算和施加，消耗了大量的前处理时间。

首先，笔者针对集中黏弹性动力人工边界前处理工作繁琐、精度不高等问题，基于一致黏弹性人工边界(CVABE)及VSBE的理论[4-5]，提出了一种全新的人工边界，即二维弧形一致黏弹性边界单元(2DACVABE)，该单元不仅具有VSBE特性，且单元刚度及阻尼矩阵与普通有限单元类似。然后，笔者深入讨论了在二维空间内2DACVABE模型中的震源频域、介质弹性波速和单元厚度对分析结果精度的影响。结果显示，2DACVABE能得到较好的精度和稳定性，并且其效果同时受到以上多个因素的影响。最后，笔者给出2DACVABE的适用范围，为进一步的研究工作提供一定的参考。

1 等效黏弹性人工边界单元

在二维一致黏弹性人工边界单元中，推导了模型的总体刚度矩阵〔式(1)〕、阻尼矩阵〔式(2)〕[5]。但是，由于总体刚度和阻尼矩阵版块在目前大型有限元分析软件设置中的缺失，模型无法模拟真实总体刚度及阻尼，因此迫切需要开发一种等效于VSBE的普通有限单元。

以往的做法是，在已建立的模型区域外表面上沿法向延伸一层同类型的单元，再固定外单元表面节点。为使边界上满足阻尼与的恢复力要求，利用矩阵等效原理(图1)，将单元的一致刚度转化为边界单元等效刚度矩阵，将阻尼矩阵转化为与边界单元等效阻尼矩阵。

图1 一致黏弹性边界单元

(1)

(2)

式中：KBT为切向刚度；KBN为法向弹簧刚度；CBT为切向阻尼系数；CBN为法向阻尼系数。

2 二维弧形一致黏弹性边界单元

2.1 二维弧形一致黏弹性边界单元的应用

等效的黏弹性边界单元回避了集中黏弹性边界设置中逐个设置弹簧-阻尼器的过程，在保证精度的前提下，简便地完成了VSBE的施加过程。但是以往的研究者在处理半无限域人工边界时，为了简化模型，多采用平直的边界。而在VSBE的推导过程中，采用的是柱面波假设，简化为平直边界，明显不符前提假设条件及实际情况，因此笔者提出一种新的黏弹性人工边界处理方式，即二维弧形一致黏弹性边界单元(2DACVABE)。

在具体设置中，首先以散射波源作用点为中心，为使得波源至弧形边界的距离不变，建立以波源为中心、半径为R的半圆形边界有限模型，然后沿边界法向在已建成的有限元模型的边界上延伸一层普通单元，固定外层边界，并通过计算赋予合适的材料系数，以替代原离散化的弹簧-阻尼器。系数的计算与设置见2.2节，模型如图2。

图2 弧形一致黏弹性边界单元

此方法无需逐一计算和布置弹簧-阻尼器，大大简化了前处理工作。同时，因模型半径为定值，避免了材料特性相关系数计算中近似取值所带来的误差，有效提高了有限元动力分析的精度。

2.2 等效刚度与等效阻尼的设置

在弧形边界单元的刚度与阻尼单元计算中，参考等效黏弹性人工边界单元的思路，使CVABE的刚度矩阵与相对应的弧形边界单元刚度矩阵相等，并采用与刚度成正比的阻尼矩阵[4]，即：

(3)

由定义可推导边界单元的等效剪切模量和等效弹性模量，如式(4)；对阻尼矩阵，为使边界单元设置更加简便，将等效单元的阻尼比例系数设置为两个方向系数的平均值[4],如式(5)：

(4)

(5)

综上所述，只需定义2DACVABE的材料特性相关系数，便可使得在有限元软件中普通单元的阻尼矩阵、等效刚度与CVABE的阻尼矩阵、刚度相等，便完成2DACVABE的设置。即边界处应力-应变情况等价于CVABE，新的方法明显简化了黏弹性边界的设置过程，并保证了较高的精度。

3 算 例

3.1 有限元模型建立

笔者在建立弧形边界单元模型时，取半径为50 m 的半圆形近场观测区域，模型周围延伸一层单位厚度单元，将外层节点固定。

选取介质密度为2 000 kg/m3、杨氏模量为20 MPa、泊松比为0.3的弹性半无限域介质进行研究。考虑介质Rayleigh阻尼，其α=0.616，β=0.000 312。2DACVABE的材性系数参考式(4)和式(5)进行计算。在中心点施加1/2个频率为1 Hz的正弦脉冲，并取对称轴处的表面点A和埋深10 m的点B为观测点进行观测，模型见图3。为验证2DACVABE的效果，同时建立固定边界、黏性边界、远置边界和集中黏弹性边界等多组模型进行对比分析。

图3 弧形一致黏弹性单元边界模型

3.2 弧形一致黏弹性边界单元精度分析

为验证2DACVABE的精度，对各边界模型观测点的位移进行对比分析(图4)，并将建模复杂、求解效率低但模拟效果最佳的远置边界作为标准解，以寻求精度接近远置边界且建模效率高、求解效率高的人工边界。

图4 观测点位移时程曲线

由图4可见：固定边界模型出现了明显的反射干扰现象，与远置边界标准解相差很大，不能反映半无限域介质中波的实际传播情况；集中黏弹性边界模型中，波的反射明显减少，但在中后段明显出现了向下漂移情况，笔者认为其精度不能满足要求；黏性边界、弧形一致黏弹性边界单元及远置边界模型在脉冲波输入后，能量迅速减小，仅产生了很小的波动，最后几乎完全被边界单元吸收，满足精度要求。将常用的人工边界与笔者提出的弧形一致黏弹性边界优劣势进行对比总结，如表1。

表1 人工边界优劣势对比

为进行人工边界的精度对比，将各边界模型观测点位移与远置边界观测点位移(标准解)相除，可得表2，并对表2进行处理，得到笔者提出的弧形一致黏弹性边界与其他常用人工边界相比后的精度提升系数，如表3。分析表2、表3可得，弧形一致黏弹性边界单元模型精度随着埋深的增加而增加。对于地表面观测点A而言，与其他常用人工边界相比，弧形一致黏弹性单元的精度提升0.6%～3.4%；对于埋深较深的观测点B而言，与其他常用人工边界相比，弧形一致黏弹性单元的精度提升3.1%～16.7%;从均值来看，弧形一致黏弹性单元的精度平均提升3.5%～15.4%。综上，笔者认为，弧形一致黏性单元的计算精度和稳定性均优于黏弹性边界和集中黏弹性边界，精度与远置边界最为接近，并且其设置方法更为简单，更易于在有限元软件中实现。

表2 时程曲线拟合精度对比

表3 弧形一致黏弹性边界单元精度提升系数

3.3 震源频率与弹性波速

首先，选取多组不同频率的动力荷载进行对比分析，以便于进一步研究震源频率与介质弹性波速对弧形一致黏弹性边界单元模型的计算精度影响。各组均为1/2个等幅值正弦脉冲波，频率分别为0.35、0.50、1.00 Hz，加载时程如图5。然后，在建模过程中，保持其他因素一致，挑选两种不同的介质弹性波速模型，波速分别为10、100 m/s。通过对频率与弹性波速的多种组合模型进行数值分析，依据图6，认为介质弹性波速和震源荷载频率均为影响2DACVABE模拟效果精确度的关键因素。

图5 脉冲波加载时程曲线

图6 观测点位移时程曲线

由图6可知，在保持介质弹性波速不变的情况下，震源频率越高，其后段反射干扰越明显，模拟效果越差；对比图6(a)、图6(c)或图6(b)、图6(d)可知，在相同波源频率的情况下，模拟效果随着介质弹性波速增大而显著改善，计算精度提高。因此，笔者得出震源的频率和介质材料性质对动力响应均有重要影响的结论。在有限元软件计算中使用2DACVABE时，建议要合理选择震源频率范围，规避频率过高的震动对结果的影响，并且要考虑材料性质对模拟效果的影响，对于介质应选取较高的弹性波速。

3.4 弧形一致黏弹性边界单元厚度

在以往的等效黏弹性边界单元研究中，研究者认为，边界单元的厚度对模拟结果影响不大，在1～50倍厚宽比范围内取值均可获得满意的结果[4]。笔者针对弧形一致黏弹性边界，选取了厚宽为0.1～5.0的5种模型进行计算，如图7。

由图7可知：当边界单元厚宽比h/l≤0.1时，出现了明显的反射干扰现象，这主要是因为单元厚宽比过小时，局部位置处可能发生边界单元的内表面触及固定外表面的情况，这就解释了为什么观测点时程图出现类似于固定边界干扰反射的情况；当h/l≥5.0时，观测点出现了明显的漂移现象，这主要是因为边界单元越厚其在相同应力下产生的法向应变越大，最终呈现为向下的漂移情况，造成了不可接受的误差。

图7 观测点位移时程曲线

综上所述，通过计算不同厚度的人工边界单元模型的曲线拟合度，发现弧形一致黏弹性边界单元厚度对动力响应有直接的影响。当1.0

4 结 论

笔者基于弧形一致黏弹性人工边界及黏弹性人工边界单元的理论，提出了一种全新的二维弧形一致黏弹性边界单元概念。具体阐述了该人工边界的建立思路和在有限元软件中的实现方法，并通过数值算例分析了其适用性：

1)考虑到集中黏弹性边界设置中近似取值所带来的误差和过于繁琐的前处理工作，笔者提出了一种新的人工边界，即二维弧形一致黏弹性边界单元。该单元模型的操作过程比集中黏弹性边界更加简单，并通过算例证明了其具有良好的精度及稳定性。

2)人工边界的动力模拟效果与输入震源的动力特性及分析介质的弹性常数存在密切的关系，建议要合理选择震源频率范围，规避频率过高的震动对结果的影响，同时要考虑材料性质对于模拟效果的影响，尽量选取较高的介质弹性波速。

3)弧形一致黏弹性边界单元厚度直接影响到模拟的计算精度水平。当厚宽比h/l超过推荐范围(1.0,2.0)时，即可能存在边界反射干扰或法向漂移等问题，使其计算结果存在较大误差甚至错误，建议合理选择单元厚度以控制误差。