地基阻抗力时域递归参数的计算方法及程序实现.pdf

第30卷 第1期 岩 土 工 程 学 报 Vol.30 No.1 2008 年 1 月 Chinese Journal of Geotechnical Engineering Jan., 2008 地基阻抗力时域递归参数的计算方法及程序实现 赵建锋 1，杜修力2 1.青岛理工大学土木工程学院，山东 青岛 266033；2.北京工业大学城市与工程安全减灾省部共建教育部重点实验室，北京 100022 摘 要基于笔者提出的考虑基础动力阻抗函数频率相关特性的地基阻抗力时域计算方法，对时域递归参数的拟合方 法和所选系统稳定条件进行了研究，编制了阻抗力时域递归参数的计算程序。基于通用有限元软件和地基阻抗力时域 差分形式，完成了考虑土–结构动力相互作用影响的结构时程分析程序二次开发。最后，通过数值算例验证了该方法 的有效性，并初步讨论了土–结构相互作用系统动力响应分析中基础动力阻抗函数频率相关性的影响。 关键词土–结构动力相互作用；动力阻抗函数；频率相关；稳定性；时域差分 中图分类号TU470 文献标识码A 文章编号1000–4548200801–0034–07 作者简介赵建锋1976– ，男，河南叶县人，讲师，博士，主要从事结构动力分析方面的研究工作。E-mail zhaojf。 Computation and realization procedure for time-domain recursive parameters of ground resistance ZHAO Jian-feng1, DU Xiu-li2 1. College of Civil Engineering, Qingdao Technological University, Qingdao 266033, China; 2. Key Lab for Urban Security foundation dynamic impedance; frequency dependent; stability; time-domain difference 0 引 言 目前我国抗震计算中常采用刚性地基假定，即把 地基或岩土看作为刚性，地震时建筑物基础的运动与 其邻近自由场地一致，未考虑土和建筑物间的相互作 用及地基辐射阻尼的影响。随着震害经验的积累和抗 震计算理论的发展，在结构动力响应分析中完整、准 确地考虑土–结构动力相互作用的影响越来越受到工 程和学术界的重视。在考虑土–结构动力相互作用影 响的结构非线性动力响应分析中，基础与结构界面间 的相互作用力的计算是关键问题之一。解决这一问题 的方法可归纳为两大类[1]①整体分析法；②子结构 法。前者可以完整、准确地考虑土体、结构介质中的 非均匀、非线性因素以及土体–结构界面间的非线性 动力相互作用，但其所需的计算成本通常是难以为工 程界接受的。后者则将结构与地基作为不同的子结构 分别进行离散和分析，然后根据交界面处的响应一致 条件加以综合，使进行非线性迭代的自由度数大为缩 减[2]。因此，子结构法受到了工程界的广泛重视。频 域子结构法能充分考虑基础阻抗函数频率相关特性， 但无法处理结构非线性问题；时–频域混合法[3]和频 –时域混合法[4-5]均需要反复付氏变换和每一计算步 的反复迭代，也加大了工作量；而时域子结构法由于 基础阻抗函数随频率变化复杂，这也就使得通过基础 阻抗函数卷积分表述的时域地基阻抗力的计算变得十 分复杂和耗时。 时域地基阻抗力的实现目前主要采用①集总参 ─────── 基金项目国家杰出青年科学基金资助项目（50325826） ；国家自然科 学基金重大研究计划项目（90715032） ；北京市教育委员会科学技术与 研究生建设项目（4700001401509） 收稿日期2006–12–04 第 1 期 赵建锋，等. 地基阻抗力时域递归参数的计算方法及程序实现 35 数模型；②脉冲反应法；③集总参数模型耦合虚拟力 法。为了简化时域地基阻抗力的计算，众多国内外学 者[6-14]已提出了各种集总参数模型。我们知道，当阻 抗函数的频率相关特性较弱时，容易将阻抗函数较为 准确地转换为时域的 Voigt 模型来表述时域地基阻抗 力。然而，在实际问题中由于土层的不均匀性以及基 础的型式和埋深等因素的影响，基础阻抗函数随频率 变化通常是很复杂的，表现出强烈的频率相关性。如 结构为线性，体系的动力分析可采用子结构法在频域 内有效地实施。但对于非线性结构，在时域上采用与 频率有关的动力阻抗进行相互作用分析就变得十分困 难与复杂。而仅采用某一频率下动力阻抗函数值的单 自由度集总参数模型又不能充分地包含所有频率成分 而难以反映实际性态。众多学者[8-14]通过增加集总参 数个数和改善组合形式来提高对阻抗函数曲线的拟合 精度。实质上，集总参数模型能很好地反映基础阻抗 函数中不具有平方可积性的奇异性特征部分，即对应 于瞬态影响的那部分。而反映延续影响、具有平方可 积性、且变化复杂的正则特征部分则未在模型中得到 完全体现，使得其计算精度的广泛适用性存在问题。 这也是集总参数模型通过引入更多参数来提高对阻抗 函数曲线的整体拟合精度效果并不明显的原因所在。 另一方面，根据傅立叶逆变换的卷积分性质和简化假 定，Wolf[15]提出一种脉冲反应法，Nakamura N[16-17] 进行了改进，其实质也是基于时频变换的积分算法， 并未根本解决时域地基阻抗力计算繁琐的问题。 Safak[18]提出将时域差分递归模型用于时域地基阻抗 力的表述，较好地模拟了基础阻抗函数中的正则项特 性，且精度在一定条件下可通过调整相应的滤波器函 数中的参数数目来控制，但由于时域差分递归模型不 能完整反映基础阻抗函数的物理实质，其相应的滤波 器函数在 Nyquist 频率附近存在本质的缺陷，这就降 低了对时域地基阻抗力的模拟精度。笔者[2]将集总参 数模型与时域差分递归的滤波器模型相结合，提出了 一种可完整考虑基础阻抗函数中奇异性特征和正则性 特征的简便、实用且精度可控的时域地基阻抗力的计 算方法，该方法将当前时刻的时域地基阻抗力的分解 表述为两项有明确物理意义的阻抗力之代数和，一是 由弹簧阻尼器集总参数模型代表的瞬态影响的阻抗力 部分，另一部分是由前几个时步的虚拟力和位移值表 示的反映延续影响的虚拟力部分。该方法在考虑土– 结构动力相互作用的结构响应分析中的应用中有两个 关键问题①选择合适方法获得拟合参数；②在动力 时程分析程序中的实现。本文针对这两个问题进行讨 论， 给出为了获得稳定的差分系统所需要的稳定条件， 然后采用最优化方法讨论了参数的有效拟合方法。最 后，在有限元软件中通过引入虚拟力的方法实现了考 虑阻抗函数频率相关性的地基阻抗力的时域计算。数 值模型结果显示，该方法简便实用，能够考虑阻抗函 数的频率相关性。 1 稳定条件 在频域中，基础–结构交界面上的相互作用力- 位移关系可写为 FSuωωω ， 1 式中， Sω为基础的动力阻抗函数， Fω和 uω分 别为相互作用力和位移幅值。通常基础动力阻抗函数 Sω可分解为奇异部分和正则部分之和，表示为 sr SSSωωω ， 2 其 中 奇 异 部 分 不 具 有 平 方 可 积 性 ， 可 写 为 s iSKCωω，正则部分是平方可积的。 根据付氏变换的线性性质，基础–结构交界面上 的相互作用力–位移关系在时域内可写为 r 0 d t F tKu tCu tS tuτττ− ∫ ， 3 其中，右端前两项相应于瞬时响应，第三项卷积相应 于延续响应。由延续响应引起的地基阻抗力可表示为 dr 0 d t FtStuF tKu tCu tτττ−−⎡⎤ ⎣⎦∫ 。 4 从式（4）可以看出如果把位移当作输入，延续 响应引起的相互作用力 d F t当作输出， 根据数字信号 处理[19]的知识，对应于阻抗函数的位移与延续响应引 起的相互作用力 d F t就构成了线性时不变（LTI）系 统（用于改变频谱形状的线性时不变系统也称为频率 成形滤波器，简称滤波器） ，而频域阻抗函数就是该 LTI滤波器系统的系统函数。 简单的离散时间滤波器可按下面的差分方程定 义 01 MN ii ii y tb x ti ta y ti t − ∆−− ∆ ∑∑ 。 5 这里tx和ty分别为原始（输入）信号和滤出 （输出）信号。 i a和 i b表示过滤因子，M和N则是滤 波器的阶数。由于仅考虑延续响应部分，式（5）中需 令 0 0b 。利用离散信号Z变换的性质，可以得到滤 波器的系统函数 1 1 1 M i i i N i i i b z Y z H z X z a z − − ⎛⎞ ⎜⎟ ⎝⎠ ∑ ∑ ， 6 式中， it ze ω∆ 。 众所周知，线性时不变系统稳定的充要条件是系 统函数的收敛域包含单位圆1z 。 又由数字信号处理 36 岩 土 工 程 学 报 2008 年 知识[19]，具有有理系统函数的因果LTI系统，只有 H z的全部极点位于单位圆内时 （也即全部极点的模 小于1）时，其收敛域包含单位圆，系统就是稳定的。 滤波器系统的极点为系统函数分母多项式方程 10 i i a z− ∑ 的根（通常以复数对的形式出现） ，滤 波器系统的稳定条件转化为系统函数分母多项式方程 的根须满足 11,2,, k zkN 。 7 在参数计算时通过在误差函数中引入惩罚项来保 证系统的极点都位于复平面单位圆内。由于采用最优 化方法来计算递归参数，选择过大的惩罚因子将影响 参数的拟合精度，而过小的惩罚因子不能保证系统的 稳定。计算结果表明在实际递归参数计算中，以未 考虑稳定性拟合时的误差函数的3～5倍作为惩罚因 子能够获得很好的计算效果。 至此，配合初始松弛条件，则由上述差分方程确 定的系统就是线性时不变和因果的。对于一个给定的 输入，即可唯一确定其输出。对于特定的基础阻抗函 数，选择一组合适的拟合参数K，C， i a和 i b使之在 整个频域的误差函数最小化 [] 2 i 1 k k VWKCHSP z ω ωωωω−≥ ∑∑ ， 8 式中， Wω为权重函数。式（8）第二项为保证系统 稳定性的附加惩罚项，P为惩罚因子， k z为滤波器系 统的极点也即方程10 i i a z− ∑ 的根。 则与之对应的 地基阻抗力就可以在时域内表达为 d F tKu tCu tF t ， 9 其中 d F t为对应延续响应的阻抗力部分 dd 11 MN ii ii F tbu ti ta F ti t − ∆−− ∆ ∑∑ ，10 可以通过虚拟力方式施加。 2 参数的计算 选择合适方法计算拟合参数是基础阻抗函数时域 化的关键问题之一，可以采用误差函数平方和的最小 化来获取参数的最佳匹配，从而将参数拟合问题转换 为函数寻优问题。为了计算简便和准确拟合参数，集 总参数模型中的刚度系数、阻尼系数与计算虚拟力的 递归参数可以统一计算。 2.1 权重函数 进行最优化参数拟合的目标函数定义为式（8） ， 其中 Wω为权重函数。权重函数的应用使得拟合函 数i KCHωω和阻抗函数 Sω在选定频率范围 的拟合精度可以进行调整。权重函数的选择需要综合 考虑输入荷载的频谱特征以及结构的自振频率，原则 上使基础位移的主频段具有更高的权重。 2.2 参数敏感性 对于任何给定的有理系统函数，有各种各样的等 效差分方程及网络结构存在，在这些众多的不同结构 中进行选择的一种重要考虑是计算的简便性。 由式 （5） 和式 （6） 所示线性时不变系统容易画出与之对应的直 接Ⅰ型信号流图，如图1所示。图1中的单位延迟系 统是由系统函数 1 z−表示的，延迟N个样本的实现采 用级联N个单位延迟。在软件实现中，N个级联的单 位延迟可按N个顺序存储寄存器予以实现。 图 1 直接Ⅰ型结构的信号流图 Fig. 1 Signal flow graph for structure of direct Ⅰ 直接Ⅰ型结构虽简单直观，但存在一些缺点。首 先是系数 i a， i b对滤波器性能的控制关系不直接，因 此调整不方便。改变任一个系数 i a将影响所有的极 点，改变任一个系数 i b将影响所有的零点。其次，对 有限字长效应太敏感容易出现不稳定现象，当阶数N 较高时其极点位置灵敏度太大。对于三阶以上的滤波 器，几乎都不采用直接型结构，而是采用其它形式的 结构。 将系统函数 H z分母多项式作因式分解， 考虑系 统函数的极点常以复数对的形式出现，再把实极点成 对组合起来，可以将系统函数 H z表示成并联形式 s 1 01 12 1 12 1 N ii i ii ee z H z a za z − −− −− ∑ 。 11 图2给出一个6阶并联型结构的系统信号流图。 并联型结构运算速度快、各基本节的误差互不影响而 且可以单独调整极点的位置，容易考虑前面所讨论的 稳定性约束条件。 事实上，并联型结构转换为直接型并不复杂且不 引入任何计算误差，就目前的计算机计算精度而言， 由直接型转换为并联型其精度也是可以保证的。考虑 到滤波器系统直接型结构对极点的敏感性，对基础阻 抗函数的参数拟合计算难以获得较为合适的结果，可 以先用滤波器的并联型拟合参数然后转换为直接型。 当然也可以采用并联型来计算虚拟力，相对于直接型 结构，计算虚拟力时采用并联型结构需要稍多的寄存 第 1 期 赵建锋，等. 地基阻抗力时域递归参数的计算方法及程序实现 37 单元和计算步。本文在编制参数拟合程序中，添加了 并联型转换为直接型的子程序。 图 2 一个复极点成对组合的 6 阶系统的并联结构 Fig. 2 Parallel- structure for a sixth-order system with the poles grouped in pairs 2.3 程序实现 应用过程中的另一个实际问题是滤波器阶数的选 择。 一般来讲滤波器的阶数越高， 拟合的精度就越好。 然而高的滤波器的阶数同时还要增加计算量，而且滤 波器的阶数增大到一定程度后，再提高对拟合精度的 影响并不大。 目标函数最小化的求解方法是一个函数寻优问 题，可以通过多种优化算法来实现，比如遗传算法、 单纯形等。本文采用遗传–单纯形[20]算法，借助 MATLAB优化工具箱提供的优化函数， 编写了计算拟 合参数K，C， i a和 i b以及拟合误差的计算程序，该 程序对目标函数增加了控制稳定性条件的惩罚项。递 归参数计算的流程如图3所示。 图 3 递归参数计算流程 Fig. 3 Calculation process of recursive parameters 3 动力时程分析程序中虚拟力的引入 由上述可以看出，基础阻抗函数的拟合参数获得 后基础阻抗力 F t的计算是完全时域的。在结构动力 时程分析中应用的关键在于每一时间步需要保存前 N步的虚拟力 d F t和位移 u t， 这需要对标准时程分 析程序做一简单修改。时域子结构法中考虑土–结构 动力相互作用影响的结构有限元动力平衡方程为 [ ] sssbsss bsbb0bbb 0 0 CCMu tu t CCcMu tu t ⎡⎤⎡⎤⎧⎫⎧⎫ ⎨⎬⎨⎬⎢⎥ ⎢⎥ ⎣⎦⎩⎭⎩⎭⎣⎦ [] sssbss bsbb0bbd 0 KKu tM x t KKku tMF t ⎡ ⎤⎧ ⎫⎧⎫⎧⎫ −− ⎨⎬⎨⎬⎨⎬⎢⎥ ⎩⎭⎩⎭⎩⎭⎣⎦ ，12 式中，下标s代表结构而b代表了刚性基础，t x 表 示地震加速度。 式中的参数分别在参数拟合结果中取 0 cC， 0 kK， d 1 M i i Ftbu ti t − ∆− ∑d 1 N i i a F ∑ ti t− ∆。 上述动力平衡方程仅比常用的结构动力时程分析 多一个虚拟力项， 完全可以采用差分法或者积分法进行 求解。作者在通用有限元程序ANSYS的基础上，完成 了能考虑土–结构动力相互作用影响的特殊求解器的 改进，使之在进行动力方程求解时自动完成虚拟力 d F t的计算与施加、动力求解、存储虚拟力 d F t和位 移 u t等步骤， 为进行土–结构系统动力响应分析中考 虑基础阻抗函数的频率相关性提供了可行的便利方法。 4 数值算例 应用本文讨论的考虑基础阻抗函数频率相关性的 地基阻抗力时域差分计算方法，对某群桩基础结构在 单向简谐荷载作用下的动力响应进行了计算和分析， 并将计算结果与采用集总参数的弹簧–阻尼模型来模 拟阻抗函数的结果进行对比分析。为了计算上的方便 对计算模型进行了一些简化处理①假定结构质量由 集中质量块表示；②土层为均匀半空间；③忽略群桩 基础转动的影响。单桩静刚度按照文献[21]提供的方 法计算获得，群桩的无量纲刚度和阻尼与无量纲频率 的关系按照文献[22，23]的方法计算。 该算例采用带刚性承台的33群桩基础的阻抗 函数， 群桩基础模型如图4（a） ， 结构模型如图4（b） ， 所采用的计算模型如图4（c） 。计算参数为土层阻 尼比 s 0.05β，土泊松比0.4ν，桩间距与桩径比 /5S d，桩土弹模比 ps /1000EE ，桩土密度比 ps /1.42ρρ，桩长细比/15L d 。 以单桩静刚度为标准进行无量纲化，该基础的水 平、竖向和摇摆阻抗函数可由表示为 [] [] [] [1] H0SHh0h0 [1] V0SVv0v0 2[1] R0SVR0R0 9i 9i i KaKkac a KaKkac a Kax Kkaca ⎫ ⎪ ⎬ ⎪ ⎭∑ ， ， ， 13 38 岩 土 工 程 学 报 2008 年 式中， [1] SH K， [1] SV K为单桩水平和竖向静刚度， 0s ad Vω 为无量纲频率， 其中 s V为土体的剪切波速，d为桩直径。 图 4 算例结构及计算模型 Fig. 4 Example system and computing model 群桩阻抗函数的实部 0 i k a和虚部 0 i c a随无量 纲频率的变化如图5中实线。采用本文讨论的集总参 数模型耦合虚拟力法，利用所编写的计算软件获得了 一组6阶参数见表1。对应的拟合结果如图5中的虚 线所示。从图中可以看出，6阶参数在较宽频段上对 频率相关阻抗函数的拟合都非常理想，事实上，更低 的阶数即可获得较为理想的拟合结果。 仅以质量块上施加竖向的频率为12.5样Hz的简 谐荷载为例，分别采用本文方法、低频平滑段拟合的 集总弹簧–阻尼模型（集总参数模型1） 、对应群桩阻 抗函数12.5 Hz处的弹簧–阻尼模型（集总参数模型 2） 计算了桩帽和质量块的位移时程， 并与以频域计算 结果通过傅里叶逆变换来获得的位移时程的精确解进 行了对比，如图6。所采用的低频平滑段集总弹簧– 阻尼模型的参数见表2，对应于群桩阻抗函数12.5 Hz 处的集总参数见表3。从图5中可以看出，无论低频 段拟合的集总弹簧–阻尼模型，还是对应群桩阻抗函 数12.5 Hz处的弹簧–阻尼模型都不能很好的模拟该 工况。只有本文方法获得的结果与精确解能够很好吻 合。采用低频段拟合获得的集总弹簧–阻尼模型无论 瞬态反应还是稳态反应均与精确解存在较大误差，尤 其是稳态反应更为严重。 而采用群桩阻抗函数12.5 Hz 处的弹簧–阻尼模型所获得的稳态反应与精确解较接 近，但是初始施加荷载时的瞬态反应存在较大误差。 只有本文方法的瞬态反应和稳态反应都与精确解很好 的吻合。表明该方法能在全频域内都能很好的拟合阻 抗函数。 图 5 无量纲化的群桩阻抗函数及拟合结果 Fig. 5 Impedance normalized pile group and parameter fitting results 表 1 六阶拟合参数 Table 1 Sixth-order fitting parameters 参数水平 竖向 摇摆 K 0.659072626 1.337524234 1.348805413 C 0.010732572 0.003216525 0.012225702 a1 1.198059803 1.523321109 0.907091301 a2 -0.01274662 0.718074419 0.481192041 a3 -0.41469954 0.029270635 -0.223223397 a4 -0.12733547 -0.20512546 -0.253722679 a5 0.019850609 -0.04850812 -0.0814001576 a6 0.004431511 0.017795586 -0.0088767072 b1 -0.68098178 -1.58850492 -0.896477439 b2 -0.20512100 -2.06295210 -0.407238898 b3 0.432842593 0.051571334 0.403596090 b4 0.027620621 0.548314362 0.164842843 b5 -0.02910144 -0.00761641 0.084010842 b6 -0.00292069 -0.02380203 0.011972951 注以 N 倍单桩静刚度进行无量纲化。 第 1 期 赵建锋，等. 地基阻抗力时域递归参数的计算方法及程序实现 39 表 2 低频平滑段集总参数模型的刚度和阻尼 Table 2 Stiffness and damping of lumped parameter model in low frequency domain 参数 刚度 阻尼 竖向 0.3225517 0.0156635 注以 N 倍单桩静刚度进行无量纲化。 图 6 简谐荷载作用下的位移时程 Fig. 6 Displacement history under simple harmonic load 表 3 对应 12.5 Hz 的集总参数 Table 3 Stiffness and damping at 12.5 Hz 参数 刚度 阻尼 竖向 0.2739727 0.026979972 注以 N 倍单桩静刚度进行无量纲化。 由上述数值结果可以看出，对于阻抗函数的实部 和虚部均接近直线的低频段，用集总参数模型来计算 地基阻抗力精度其精度足够工程应用。但考虑对于诸 如近断层地震动、爆炸、冲击、动力机械装置等包含 丰富高频成分的荷载时，需要考虑阻抗函数的频率相 关性。在进行结构动力响应分析时，如果考虑结构的 非线性因素和土–结构动力相互作用的影响，通常的 时域方法无法处理地基阻抗函数的频率相关性，而频 域方法难以结构的非线性因素，都需要进行线性等效 化处理。本文提供了一种方便的可用于时域考虑地基 阻抗函数频率相关性的地基阻抗力计算方法，为土– 结构相互作用系统的动力响应准确分析提供了可能。 5 结 语 笔者[2]将集总参数模型与时域差分递归的滤波器 模型相结合，提出了一种可完整考虑基础阻抗函数中 奇异性特征和正则性特征的简便、实用的时域地基阻 抗力的计算方法。该方法在考虑土–结构动力相互作 用的结构响应分析中应用的两个关键问题是地基阻抗 力时域递归参数的计算和该方法在动力时程分析程序 中的实现。本文对此进行分析，给出为了获得稳定的 差分系统所需要的条件，即系统函数的所有极点位于 复平面单位圆内 k 1 1,2,,zkN。然后讨论了 参数的最优化拟合方法，即首先采用并联结构的系统 函数来拟合基础阻抗函数， 然后转换为等效的直接型。 拟合过程以误差函数为目标函数，通过引入惩罚项保 证了该方法的稳定性，采用遗传-单纯形算法，借助 MATLAB优化工具箱提供的优化函数， 编写了参数拟 合的计算程序。 数值算例结果验证了该方法的有效性、 简便性，能够在全频域内很好的考虑阻抗函数的频率 相关性。在此基础上初步讨论了土–结构相互作用系 统动力响应分析中基础动力阻抗函数频率相关性对结 构动力响应的影响。 参考文献 [1] WOLF J P. Dynamic soil-structure interaction[M]. New Jersey Prentice-Hall, 1985. [2] 赵建锋, 杜修力. 精度可控地基阻抗力的一种时域计算方 法[C]//第七届全国土动力学学术会议. 北京 清华大学出 版社, 2006 567–577. ZHAO Jian-feng, DU Xiu-li. An accuracy controllable time-domain approach to calculate foundation resisting force[C]// The 7th National Conference on Soil Dynamics. Beijing Tsinghua University Press, 2006 567–577. in Chinese [3] BERNAL D, YOUSSEF A. A hybrid time frequency domain ulation for non-linear soil structure interaction[J]. Earthquake Engineering and Structural Dynamics, 1998, 277 673–685. [4] DARBRE G R. Seismic analysis of non-linearly base-isolated soil-structure interacting reactor building by way of the hybrid-frequency-time-domain procedure[J]. Earthquake Engineering and Structural Dynamics, 1990, 19 725–738. [5] DARBRE G R. Application of the hybrid frequency–time- domain procedure to the soil–structure interaction analysis of a shear building with multiple nonlinearities[C]// Proceedings of the 10th World Conference Earthquake Engineering, Madrid, Spain. 1992 442–454. [6] LYSMER J, RICHART F E. Dynamic response of footings to vertical loading[J]. J Soil Mech Found Div, ASCE, 1966, 921 65–91 [7] 金 峰, 张楚汉, 王光纶. 拱坝–地基动力相互作用的时 40 岩 土 工 程 学 报 2008 年 域模型[J]. 土木工程学报, 1997, 301 43–50. JIN Feng, ZHANG Chu-han, WANG Guang-lun. A time domain model of arch dam-rock foundation interaction[J]. China Civil Engineering Journal, 1997, 301 43–50. in Chinese [8] WOLF J P, SOMAINI D R. Approximate dynamic model of embedded foundation in time domain[J]. Earthquake Engineering and Structural Dynamics, 1986, 145 683–703 [9] BARROS DE, LUCO J E. Discrete models for vertical vibrations of surface and embedded foundation[J]. Earthquake Engineering and Structural Dynamics, 1990, 192 289–303. [10] JEAN W Y, LING T W, PENZIEN J. System parameter of soil foundation for time domain dynamic analysis[J]. Earthquake Engineering and Structural Dynamics, 1990, 19 541–553. [11] WU W H, LEE W H. Systematic lumped-parameter models for foundations based on polynomial-fraction approximation[J]. Earthquake Engineering and Structural Dynamics, 2002, 317 1383–1412. [12] WOLF J P. Spring-dashpot-mass models for foundation vibrations[J]. Earthquake Engineering and Structural Dynamics, 1997, 26 931–949. [13] 栾茂田, 林 皋. 地基动力阻抗的双自由度集总参数模型 [J]. 大连理工大学学报, 1996, 364 477–482. LUAN Mao-tian, LIN Gao. 2-DOF lumped-parameter model of dynamic impedances of foundation soils[J]. J of Dalian University of Technology, 1996, 364 477 – 482. in Chinese [14] 侯兴民, 廖振鹏. 表面矩形基础阻抗函数的集中参数模型 [J]. 地震工程与工程振动, 1999, 194 6–13. HOU Xing-min, LIAO Zhen-peng. Lumped parameter models for impedance matrix of rigid rectangular foundations of the half-space[J]. Earthquake Engineering and Engineering Vibration, 1999, 194 6–13. in Chinese [15] WOLF J P, OBERNHUBER P. Nonlinear soil-structure interaction analysis using dynamic stiffness or flexibility of soil in the time domain[J]. Earthquake Engineering and Structural Dynamics 1985, 13 195 –212. [16] NAKAMURA N. A practical to trans frequency dependent impedance to time domain[J]. Earthquake Engineering and Structural Dynamics, 2006, 352 217– 231. [17] NAKAMURA N. Improved s to trans frequency-dependent complex stiffness to time domain[J]. Earthquake Engineering and Structural Dynamics, 2006, 358 1037–1050. [18] SAFAK E. Time-domain representation of frequency- dependent foundation impedance functions[J]. Soil Dynamics and Earthquake Engineering. 2006, 261 65–70. [19] OPPENHEIM A V, SCHAFER R W, BUCH J R. Discrete-time signal processing[M]. 2nd ed. New Jersey Prentice-Hall Inc, 1999. [20] 韩 炜, 廖振鹏. 一种全局优化算法遗传算法–单纯形 法[J]. 地震工程与工程振动, 2001, 212 6–12. HAN Wei, LIAO Zhen-peng. A global optimization algorithm genetic algorithm-simplex[J]. Earthquake Engineering and Engineering Vibration, 2001, 212 6–12. in Chinese [21] NOVAK M, SHARNOUBY B E L. Stiffness constants of single piles[J]. Journal of Geotechnical Engineering, ASCE, 1983, 1097 961–974. [22] MAKRIS N, GAZETAS G. Dynamic pile-soil-pile interaction Part II Lateral and seismic response[J]. Earthquake Engineering and Structural Dynamics, 1992, 212 145– 162. [23] MAKRIS N, GAZETAS G. Displacement phase differences i