The Application of Differential Quadrature Method in Structural Dynamic Analysis
-
摘要: 微分求积法（DQM）是1种求解微分方程初（边）值问题的数值方法，通常以较小的计算工作量即可获得较高的数值精度。这种方法应用于工程领域时多用来解决梁、板等结构的静力分析或结构特征值分析等问题，即对边值问题的微分方程的求解。结构动力分析属于初值问题，荷载和结构反应都具有特殊性，直接套用DQM求解边值问题并不能获得问题的解。本文尝试利用微分求积原理建立求解结构动力反应的具体方法。借鉴单元法的思想，将荷载持时划分为若干个时步，在每个时步内对动态荷载和结构反应进行离散，然后用DQM对时步逐个进行求解，得到体系在整个时域内的反应过程。通过对3种不同自振周期的线弹性单自由度体系在不同频率简谐激励下反应的计算，阐释了本文方法的可行性以及高精度、高效率的特点，通过数值试验确定了时步内相对较优的节点数，并为时步长度的选取提供了建议。Abstract: The differential quadrature method (DQM) is a numerical technique of solving initial/boundary value problems of differential equations and capable of obtaining a higher numerical accuracy with a smaller calculation workload. This method is often used to solve the problems of structural static analysis of beams and slabs or eigenvalue analysis when it is applied to the engineering fields, which is to solve the differential equation of the boundary value problems. Dynamic analysis of structures is an initial value problem as well as particular loads and structural response. As a result, applying the DQ method of solving the boundary value problem directly cannot obtain solution of problem. The principle of differential quadrature is applied to establish the specific method of solving structural dynamic response in this paper. By analogizing the idea of unit method, the duration of the load is divided into many time steps and dynamic load and structure are discretized in each time step, then the response of the system can be solved in the whole time domain by employing the DQ method step by step. The feasibility of this method and the characteristics of high precision and high efficiency are expounded by calculating the response of three linear elastic single-freedom-degree system of different natural vibration periods excited by simple harmonic loads of different frequencies. By means of numerical experiment, the optimal meshing scheme is determined and the suggestion for the time step is given.
-
引言
地震等灾害环境的作用使结构产生复杂的动态响应，可能导致结构失效、破坏甚至倒塌，从而形成灾害、造成损失。对新建结构进行抗灾设计及对已建结构进行抗灾加固是防御和减轻工程灾害及其损失的有效途径，这就要求认识结构在这些灾害作用下的性态和响应行为，而结构动力分析则是实现这一目标的基本手段。
结构动力分析的本质是实现对动力荷载激励下的结构运动微分方程（组）的求解，属于常微分方程的初值问题。目前用于分析结构动力的方法大致可以分为2类：一类是变换方法，如振型叠加法利用振型的正交性和完备性将结构动态响应向各阶振型分解，再通过叠加各阶振型响应以获得结构动态响应的结果，这类方法采用叠加原理，一般只用于线弹性结构动力分析，且要求结构具有经典阻尼特性；另一类为直接方法，即直接对结构运动微分方程进行求解而不必引入任何假定，这类方法既可用于线弹性结构，也可用于非线性结构的分析，以Newmark-β法（ Newmark，1959）等逐步积分法为代表，一般采用数值求解手段。
由于现代结构不断向大型化、复杂化发展，加之结构精细化模型的采用，导致结构动力分析的计算需求呈爆发式增长，因而需要寻求高效率的分析方法。微分求积法（Differential Quadrature Method，DQM）是由Bellman等（1971，1972）发展起来的1种求解微分方程的数值方法，可通过较小的计算工作量获得较高的计算精度。DQM要求求解域比较规则，在工程分析中一般用于时间无关问题的求解。如在结构力学领域，多用来求解静力问题和固有振动问题等（Bert等，1988, 1993, 1996, 1997；Jang等，1989；Kang等，1995, 1996；Liew等，1996a, 1996b；Sherbourne等，1991；Striz等，1988；Wang等，1993, 1994；Zeng等，2001）。这类问题都属于边值问题，即在DQM中计算的是解函数对空间坐标的导数。而结构动力分析属于初值问题，使用DQM对其求解的研究工作相对较少，Fung（2001a，2001b）、Liu等（2008）、李鸿晶等(2011a, 2011b)和廖旭等（2013）开展过相关的研究工作。本文在此基础上发展1种基于DQM的结构动力分析的高精度方法。不同于静力边值问题，实际结构的动力响应问题有其特殊性，许多情况下时间跨度很长，像边值问题一样一次性对所有时间区域进行离散求解，将出现病态问题而产生错误。本文借鉴单元法的思想，以期提高结构动力分析的计算效率。
1. 微分求积法基本原理
微分求积法是1种用于求解微分方程的数值方法。它的实质是将函数在某一离散节点处的各阶导数值，近似表示成计算域内所有节点处离散函数值的线性加权和，从而将复杂的微分方程化为关于离散点的线性方程（组）。由于本文讨论的是结构动力反应的计算，求解的运动微分方程仅是关于时间t的常微分方程，因此仅介绍一维区域内的微分求积原理。
设函数f(x)为在区间[a，b]上k阶连续可微，将区间[a，b]划分为m段，共（m+1）个互不相同的节点，分别记为x0，x1，……，xm-1，xm，其中x0=a，xm=b。
根据计算数学的函数逼近理论，函数f（x）可做如下逼近：
$$ f(x) \approx \sum\limits_{j = 0}^m {{q_j}(x)f({x_j})} $$ (1) 其中，$ {q_j}(x)$为函数空间中各线性无关的基函数。
对式（1）求k阶导数，然后将所有节点代入，得到：
$$ {f^{(k)}}({x_i}) \approx \sum\limits_{j = 0}^m {q_j^{(k)}({x_i})f({x_j}){\rm{ }}}, \;\;i = 0, \;1, \cdots \cdots, m $$ (2) 记$a_{ij}^{(k)} = q_j^{(k)}({x_i}), {f^{(k)}}({x_i}) = f_i^{(k)}, f({x_j}) = {f_j} $，则：
$$ f_i^{(k)} \approx \sum\limits_{j = 1}^N {{a^{(k)}}_{ij}{f_j}{\rm{ }}}, \;\;i = 0, \;1, \cdots \cdots, m $$ (3) 式（3）即为一维区域微分求积的基本公式。
微分求积法的关键问题是确定式（3）中的系数$ {a^{(k)}}_{ij}$。该系数称为微分求积权系数，它与网格节点的分布及选择的基函数空间类型有关，当节点的位置确定后，将基函数代入式（3），即可得到相应的权系数。
常用的函数空间是（m+1）维多项式空间，选择该空间中的所有基函数都能得到相同的权系数。最常见的多项式基函数是幂指数插值函数和拉格朗日插值函数2种，分别如式（4）、（5）所示：
$$ f(x) = {x^r}, r = 0, 1, \cdots \cdots, m $$ (4) $$ f(x) = \prod\limits_{i = 0\atop i \ne j}^m {\frac{{x - {x_i}}}{{{x_j} - {x_i}}}} \;\;{\rm{, }}\;\;j = 0, 1, \cdots \cdots, m $$ (5) 选用式（4）的幂指数函数计算权系数，得到权系数的隐式表达式，需要求解范德蒙矩阵的逆矩阵，但当m较大时，不但计算量大，而且矩阵将容易出现病态。为了解决这些问题，目前大多采用式（5）给出的拉格朗日插值基函数，因为拉格朗日插值的各项系数为各节点的函数值，形式与式（1）完全相同，直接求k阶导数，便可得到相应的k阶权系数的显示表达式，计算效率大幅度提高。
实际计算权系数时往往只需计算1阶导数的权系数，其它各阶导数的权系数可由1阶权系数以矩阵的形式方便地表示：
$$ {{\bf{A}}^{k}}={({\bf{A}})^k} $$ (6) 其中，A表示1阶权系数矩阵，A（k）表示k阶权系数矩阵。$ {\bf{A}} = \left[ {\begin{array}{*{20}{c}} {{a_{00}}}&{{a_{01}}}& \cdots &{{a_{0m}}}\\ {{a_{10}}}&{{a_{11}}}& \cdots &{{a_{1m}}}\\ \vdots&\vdots&\ddots&\vdots \\ {{a_{m0}}}&{{a_{m1}}}& \cdots &{{a_{mm}}} \end{array}} \right]$，$ {{\bf{A}}^{\left(k \right)}} = \left[ {\begin{array}{*{20}{c}} {{a^{(k)}}_{00}}&{{a^{(k)}}_{01}}& \cdots &{{a^{(k)}}_{0m}}\\ {{a^{(k)}}_{10}}&{{a^{(k)}}_{11}}& \cdots &{{a^{(k)}}_{1m}}\\ \vdots&\vdots&\ddots&\vdots \\ {{a^{(k)}}_{m0}}&{{a^{(k)}}_{m1}}& \cdots &{{a^{(k)}}_{mm}} \end{array}} \right]$, 1阶权系数${a_{ij}} = {a^{\left( 1 \right)}}_{ij} $，可通过选定拉格朗日基函数，直接得到如下的表达式：
$$ {a_{ij}} = \left\{ \begin{array}{l} \frac{{M({x_i})}}{{({x_i} - {x_j})M({x_j})}}\;\;\;, \;\;\;i \ne j\\ - \sum\limits_{k = 0\atop k \ne i}^m {{a_{ik}}} \;\;\;, \;\;\;i = j \end{array} \right. $$ (7) 其中，$ M\left({{x_i}} \right){\rm{ = }}\prod\limits_{k = 0\atop k \ne i}^m {({x_i} - {x_k})} {\rm{ }}$，$M({x_j})= = \prod\limits_{k = 0\atop k \ne j}^m {({x_j} - {x_k})} $。
应用微分求积法时还需确定采样网格节点的位置，大致分为均匀网格点和非均匀网格点2大类。虽然均匀网格点的精度总体上没有非均匀点高，但是其在处理离散荷载，如地震荷载或风荷载时，可直接使用原始采样点作为节点，不需要对荷载进行额外的插值，计算效率比不均匀网格点更高。
2. 结构动力反应微分求积分析方法
用上述微分求积法求解结构的动力反应，并以线弹性单自由度体系为例进行讨论。虽然实际结构大都为多自由度体系，且为非线性体系，但其动力反应都可以通过线性迭代和振型分解转化为线弹性单自由度体系。因而，讨论线弹性单自由度体系的微分求积法更具普遍意义，且简单直观、易于理解。
线弹性单自由度体系的动力反应运动方程为：
$$ \frac{{{{\rm{d}}^2}u}}{{{\rm{d}}{t^2}}} + 2\xi \omega \frac{{{\rm{d}}u}}{{{\rm{d}}t}} + {\omega ^2}u = p(t) $$ (8) 其中，ω、ξ分别表示体系的自振频率和阻尼比；u表示体系的位移；$ \frac{{{\rm{d}}u}}{{dt}}$、$ \frac{{{{\rm{d}}^2}u}}{{d{t^2}}}$分别表示体系的速度和加速度；p（t）为结构所受到的随时间变化的荷载。
在实际工程中，动荷载随类型的不同，作用时间差异很大，短则瞬间（如冲击荷载），长则几分钟甚至几十分钟（如风荷载）。对于一些长时间作用的荷载，在整个荷载作用时域内实施微分求积法，要保证计算结果的精确性，需要成千上万的时间节点，这将导致求解时系数矩阵的阶数过于庞大，引起矩阵的条件数大，病态效果严重，无法得到可靠的结果。为了解决此问题，可以借鉴有限单元法的思想，将整个荷载持时划分为多个等间距的时间单元，称为1个时步，在每个时步内，使用微分求积法求解。
考虑动力荷载持时内长度$ \Delta t$的时步[tj，tk]内的反应，tj、tk一般与荷载p（t）的采样时刻点重合。在时步内定义一局部坐标$tt \in {\rm{[0, }}\;\Delta t{\rm{]}} $，坐标起点为tj，方向同时间坐标。为了方便处理，正则化局部坐标，作$\tau= tt/\Delta t $，则定义域被正则化为[0, 1]，时步内关于$\tau $的运动方程可写为：
$$ \ddot u(\tau) + 2\xi \omega \Delta t\dot u(\tau) + {(\Delta t)^2}{\omega ^2}u(\tau) = {(\Delta t)^2}p(\tau \Delta t) $$ (9) 其中，$\dot u(\tau) $和$\ddot u(\tau)$分别表示$ u(\tau)$对局部坐标τ的1阶和2阶导数。
将时步离散为m段，记节点为τ1，τ2，……，τm，将$ u(\tau= {\tau _i})$，$ \dot u(\tau = {\tau _i})$，$ \ddot u(\tau = {\tau _i})$，$p(\tau = {\tau _i}) $分别简记为ui，${\dot u_i} $，$ {\ddot u_i}$，pi，则在各离散节点处的方程为：
$$ {\ddot u_i} + 2\xi \omega \Delta t{\dot u_i} + {(\Delta t)^2}{\omega ^2}u{}_i = {(\Delta t)^2}p{}_i, \;\;i = 0, \;1, \cdots \cdots, m $$ (10) 对这些节点使用微分求积法：
$$ {\dot u_i} = \sum\limits_{j = 0}^N {{a_{ij}}{u_j}}, \;\;i = 0, \;\;1, \cdots \cdots, \;m $$ (11) $$ {\ddot u_i} = \sum\limits_{k = 0}^m {{a_{ik}}{{\dot u}_k} = } \sum\limits_{k = 0}^m {{a_{ik}}\sum\limits_{j = 0}^m {{a_{kj}}{u_j}} }, \;\;i = 0, \;\;1, \cdots \cdots, \;m $$ (12) 其中，${a_{ij}}$是微分求积的权系数，将式（11）、（12）写成矩阵的形式：
$$ {\bf{\dot u}}={\bf{Au}} $$ (13) $$ {\bf{\ddot u}} = {{\bf{A}}^2}{\bf{u}} $$ (14) 其中，$ {\bf{u}} = {({u_0}, {u_1}, \ldots \ldots, {u_m})^{\rm{T}}}$。
已知时步的初始位移和速度分别为u0和v0，且$ {\dot u_0} = {v_0}\Delta t$，将式（11）、（12）改写为：
$$ {{\bf{\dot u}}_{\bf{s}}} = {{\bf{G}}_{\bf{0}}}u{}_0 + {\bf{G}}{{\bf{u}}_{\bf{s}}} $$ (15) $$ {{\bf{\ddot u}}_{\bf{s}}} = {{\bf{G}}_{\bf{0}}}v{}_0\Delta t + {{\bf{G}}_{\bf{0}}}{\bf{G}}{u_0} + {{\bf{G}}^{\bf{2}}}{{\bf{u}}_{\bf{s}}} $$ (16) 其中，$ {{\bf{u}}_{\bf{s}}} = {({u_1}, {u_2}, \ldots \ldots, {u_m})^{\rm{T}}}$，$ {{\bf{G}}_{\bf{0}}} = ({a_{10}}, {a_{{\rm{2}}0}}, \ldots \ldots, {a_{m0}})$，$ {\bf{G}} = \left[ {\begin{array}{*{20}{c}} {{a_{11}}}& \cdots &{{a_{m1}}}\\ \vdots&\ddots&\vdots \\ {{a_{1n}}}& \cdots &{{a_{mm}}} \end{array}} \right]$。
将式（15）、（16）代入式（10），并将已知量与未知量分离在等式的两侧，整理后得到如下方程：
$$ {\bf{\hat k}}{{\bf{u}}_{\bf{s}}}={\bf{\hat p}} $$ (17) 其中，${\bf{\hat k}} = {{\bf{G}}^{\bf{2}}} + 2\xi \omega \Delta t{\bf{G}} + {(\Delta t)^2}{\omega ^2}{\bf{I}} $，${\bf{\hat p}} = {(\Delta t)^2}{\bf{p}} - {u_0}(2\xi \omega \Delta t + {{\bf{G}}_{\bf{s}}}{{\bf{G}}_{\bf{0}}}) - {v_0}\Delta t{{\bf{G}}_{\bf{0}}} $，I是m阶单位矩阵，$ {\bf{p}} = {({p_1}, {p_2}, \ldots \ldots, {p_m})^{\rm{T}}}$。
式（17）为各时步内微分求积的基本方程，求解该方程可得到时步内各点的位移反应，再运用微分求积原理，可进一步求得各节点速度反应：
$$ {{\bf{v}}_{\bf{s}}} = \frac{1}{{\Delta t}}({{\bf{G}}_{\bf{0}}}{u_0} + {\bf{G}}{{\bf{u}}_{\bf{s}}}) $$ (18) 其中，${{\bf{v}}_{\bf{s}}} = {\left({{{\left. {\frac{{{\rm{d}}u}}{{{\rm{d}}t}}} \right|}_{\tau = {\tau _1}}}, {{\left. {\frac{{{\rm{d}}u}}{{{\rm{d}}t}}} \right|}_{\tau = {\tau _2}}}, \ldots \ldots, {{\left. {\frac{{{\rm{d}}u}}{{{\rm{d}}t}}} \right|}_{\tau = {\tau _m}}}} \right)^{\rm{T}}} $。
将本时步末的位移$ {u_m}$和速度${\left. {\frac{{{\rm{d}}u}}{{{\rm{d}}t}}} \right|_{\tau = {\tau _m}}} $作为下一时步的初始位移和速度，从初始0时刻开始逐时步进行求解，可得到荷载作用的所有时刻的位移和速度反应。除了位移和速度外，工程上关心的加速度$ \frac{{{{\rm{d}}^2}u}}{{{\rm{d}}{t^2}}}$可通过运动方程（8）变形得到：
$$ \frac{{{{\rm{d}}^2}u}}{{{\rm{d}}{t^2}}} = p(t) - 2\xi \omega \frac{{{\rm{d}}u}}{{{\rm{d}}t}} - {\omega ^2}u $$ (19) 3. 算例
通过具体的算例验证上述微分求积法求解结构动力反应计算结果的精确性和可靠性。选择3种不同自振周期的单自由度体系，其频率范围大致覆盖低频、中频和高频，阻尼比都取为工程中常见的0.05，在上面施加不同频率的正弦荷载，体系的基本参数如表 1所示，荷载信息如表 2所示。
表 1 体系的基本特性Table 1. Basic characteristics of systems体系编号 自振周期Tn/s 阻尼比ξ 1 0.68 0.05 2 0.25 0.05 3 0.08 0.05 表 2 简谐荷载的信息Table 2. Information of simple harmonic load荷载函数/N 周期/s 振幅/m 持时/s sin2πt 1 1 40 sin10πt 0.2 1 40 sin20πt 0.1 1 40 用微分求积法求解3种体系在以上简谐荷载下的动力反应。为了方便计算，时步内采用均匀网格离散方案，时步的长度Δt取为简谐荷载的周期。由于1个周期的长度是简谐荷载的最小重复单元，这样取值不仅能方便地分析时步分段数m对数值稳定性和精度的影响，更能清楚地发现荷载周期与步长Δt取值的规律。
首先，计算该时步长度下采用不同m所得到的体系的位移和速度反应，然后采用体系在简谐荷载下的解析解作为精确解进行校核。表 3—5列出了分段数m为20内的偶数时，简谐荷载激励下3种体系动力反应DQM解与精确解的平均相对误差。
表 3 荷载周期1s的平均相对误差Table 3. Average relative error with load period of 1sm 体系编号 1 2 3 位移/% 速度/% 位移/% 速度/% 位移/% 速度/% 2 98.734 98.734 98.734 98.734 98.734 98.734 4 88.218 81.849 99.415 389.850 53.532 1293.500 6 42.402 32.547 15.926 4.728 10.979 9.609 8 6.698 8.020 1.543 7.814 2.042 2.741 10 1.623 2.687 0.907 2.338 1.036 0.349 12 0.121 0.138 0.897 3.341 0.819 4.156 14 0.008 0.009 0.935 5.153 0.716 1.290 16 0 0 0.930 3.899 0.635 0.825 18 0 0 1.005 3.292 0.577 2.28 20 0 0 1.268 4.715 0.527 0.899 表 4 荷载周期0.2s的平均相对误差Table 4. Average relative error with load period of 0.2sm 体系编号 1 2 3 位移/% 速度/% 位移/% 速度/% 位移/% 速度/% 2 99.749 99.749 99.749 99.749 99.749 99.749 4 799.760 2460.400 403.870 411.780 108.540 225.480 6 90.884 11.229 23.620 19.707 11.436 5.439 8 6.052 15.762 2.102 2.135 1.068 2.520 10 0.264 0.046 0.116 0.140 0.765 1.155 12 0.009 0.022 0.005 0.004 0.807 3.325 14 3681.600 256.240 251.760 2935.800 972.740 1334.400 16 0 0 0 0 0.093 0.167 18 0 0 0 0 0.010 0.013 20 0 0 0 0 0.001 0.001 表 5 荷载周期为0.1s的平均相对误差Table 5. Average relative error with load period of 0.1 secondsm 体系编号 1 2 3 位移/% 速度/% 位移/% 速度/% 位移/% 速度/% 2 99.875 99.875 99.875 99.875 99.875 99.875 4 1620.200 1370.200 578.690 459.540 51.868 52.298 6 175.570 9.444 65.345 7.903 6.80×10282 8.63×10282 8 11.571 9.569 4.413 3.254 1.951 2.655 10 0.498 0.046 0.204 0.044 0.142 0.266 12 0.016 0.013 0.007 0.005 0.006 0.007 14 27683 306.960 1322.400 246.430 234.600 3888.600 16 0 0 0 0 0 0 18 0 0 0 0 0 0 20 0 0 0 0 0 0 从表 3—5可以看出，当时步分段数m很小时，DQM计算结果严重偏离精确解，但随着m的增大，除少部分点外，误差大致呈迅速减小的趋势，个别数据（如m=14时）反应误差巨大是由于算法在该时步分段下，时步长度取为荷载周期时出现了数值不稳定现象，初始误差随着逐时步推进不断放大，最后完全湮灭真实的结果。当m增大到10时，3种体系在不同周期的简谐荷载下的位移和速度反应都减小到了5%以内，对于实际工程已足够精确。为了更形象地描述这种精确程度，图 1—3给出了0—20s内荷载周期为1s下，m=10时DQM的位移计算结果与精确解，为了更清晰地进行对比，对每张图进行了局部的放大。
图 1 简谐荷载激励下体系1的位移反应Figure 1. Displacement response of system 1 under simple harmonic load
图 2 简谐荷载激励下体系2的位移反应Figure 2. Displacement response of system 2 under simple harmonic load
图 3 简谐荷载激励下体系3的位移反应Figure 3. Displacement response of system 3 under simple harmonic load从图中可以看出，DQM的计算结果与精确解几乎完全重合，充分说明了分段数m=10时，DQM计算动力反应足够精确，且数值稳定性能得到很好的保证。以地震荷载为例，中低频地震荷载的大部分周期一般在0.2s以上，偏于保守取为0.2s，分10段后节点间距为0.02s。而目前一般的地震地面加速度采样周期仅为0.005s或0.1s，节点间距取采样周期的几倍以上，都能得到精确的结果，计算精度较高。
对比表 3、4、5可以发现，当m由10继续向上增大时（排除m=14时因失稳误差剧增的情况），误差基本上呈继续减小的趋势，很多情况下误差甚至降至1%或1‰以内。但这对工程实际已无太大的意义，反而会随着m的增大，不断地增大计算量。因此，对于结构在简谐荷载作用下的反应，采用均匀节点方案，即在1个荷载周期内取m=10，时步内插入9个内节点，既能达到土木工程领域内精度的要求，又能最大限度地减少计算工作量，是相对较优的时步节点数量。而且，以往DQM求解工程结构静力问题的大量经验表明，m取8—10可得到满意的计算结果（王鑫伟，1995），这与本文得出的m=10的取值较优也相符，这虽然是DQM解决静力问题的经验，但在1个周期内的动力问题与静力问题存在不少的联系，从而也可间接说明本文将m=10作为较优参数的合理性。
由表 3—5已在工程可接受的范还可以发现，对于简谐荷载激励下的反应，当m=10时，取时步长度Δt与荷载周期相等时，计算误差已在工程可接受的范围内。因此，计算体系在简谐荷载激励下的反应，在选定m=10的前提下，将时步长度和简谐荷载的周期保持一致即可。
虽然简谐荷载是1种理想的荷载，但是对计算实际荷载激励下的反应具有重要的意义，这是因为工程中实际的动荷载都是一定持时的暂态荷载，都可以通过周期延拓表示成一系列简谐荷载的叠加。实际计算时，可以首先采用谐波分析，检测出其包含的不同周期的简谐波的幅值，得到该动力荷载的频谱；然后以最大幅值所对应的周期，即卓越周期为基准，确定等效周期；最后，取时步的长度为荷载的等效周期，将时步等分成10段进行计算，即可得到较精确的计算结果。
4. 结论
本文将微分求积法引入结构动力反应的分析与计算，并通过数值算例得到如下的结论：
（1）采用微分求积法求解结构在动力荷载激励下的反应合理可行。在较大的节点距离下依然能够得到精确的结果，计算精度高，且具有普适性，对不同自振周期的结构、不同频率的动力荷载都适用。
（2）用DQM进行动力分析时，能一次性求得多个时刻的反应。相比传统的每次只能求1个时刻的逐步积分法，计算效率得到提高，计算成本也得到了降低。
（3）在时步长度Δt一定时，DQM求解动力反应的计算精度和数值稳定性与时步分段数m有关。排除一些失稳飘移的情况，一般m越大，计算精度越高，但计算量也越大。综合考虑，对于均匀网格离散方案，实际计算时取m=10相对较优。
（4）使用DQM进行实际结构动力反应分析时，时间步长Δt可选为动力荷载的等效周期，然后将各时步等分成10段来计算，这样可获得较满意的计算结果。
-
图 1 简谐荷载激励下体系1的位移反应
Figure 1. Displacement response of system 1 under simple harmonic load
图 2 简谐荷载激励下体系2的位移反应
Figure 2. Displacement response of system 2 under simple harmonic load
图 3 简谐荷载激励下体系3的位移反应
Figure 3. Displacement response of system 3 under simple harmonic load
表 1 体系的基本特性
Table 1. Basic characteristics of systems
体系编号 自振周期Tn/s 阻尼比ξ 1 0.68 0.05 2 0.25 0.05 3 0.08 0.05 
下载: 导出CSV
表 2 简谐荷载的信息
Table 2. Information of simple harmonic load
荷载函数/N 周期/s 振幅/m 持时/s sin2πt 1 1 40 sin10πt 0.2 1 40 sin20πt 0.1 1 40
下载: 导出CSV
表 3 荷载周期1s的平均相对误差
Table 3. Average relative error with load period of 1s
m 体系编号 1 2 3 位移/% 速度/% 位移/% 速度/% 位移/% 速度/% 2 98.734 98.734 98.734 98.734 98.734 98.734 4 88.218 81.849 99.415 389.850 53.532 1293.500 6 42.402 32.547 15.926 4.728 10.979 9.609 8 6.698 8.020 1.543 7.814 2.042 2.741 10 1.623 2.687 0.907 2.338 1.036 0.349 12 0.121 0.138 0.897 3.341 0.819 4.156 14 0.008 0.009 0.935 5.153 0.716 1.290 16 0 0 0.930 3.899 0.635 0.825 18 0 0 1.005 3.292 0.577 2.28 20 0 0 1.268 4.715 0.527 0.899
下载: 导出CSV
表 4 荷载周期0.2s的平均相对误差
Table 4. Average relative error with load period of 0.2s
m 体系编号 1 2 3 位移/% 速度/% 位移/% 速度/% 位移/% 速度/% 2 99.749 99.749 99.749 99.749 99.749 99.749 4 799.760 2460.400 403.870 411.780 108.540 225.480 6 90.884 11.229 23.620 19.707 11.436 5.439 8 6.052 15.762 2.102 2.135 1.068 2.520 10 0.264 0.046 0.116 0.140 0.765 1.155 12 0.009 0.022 0.005 0.004 0.807 3.325 14 3681.600 256.240 251.760 2935.800 972.740 1334.400 16 0 0 0 0 0.093 0.167 18 0 0 0 0 0.010 0.013 20 0 0 0 0 0.001 0.001
下载: 导出CSV
表 5 荷载周期为0.1s的平均相对误差
Table 5. Average relative error with load period of 0.1 seconds
m 体系编号 1 2 3 位移/% 速度/% 位移/% 速度/% 位移/% 速度/% 2 99.875 99.875 99.875 99.875 99.875 99.875 4 1620.200 1370.200 578.690 459.540 51.868 52.298 6 175.570 9.444 65.345 7.903 6.80×10282 8.63×10282 8 11.571 9.569 4.413 3.254 1.951 2.655 10 0.498 0.046 0.204 0.044 0.142 0.266 12 0.016 0.013 0.007 0.005 0.006 0.007 14 27683 306.960 1322.400 246.430 234.600 3888.600 16 0 0 0 0 0 0 18 0 0 0 0 0 0 20 0 0 0 0 0 0
下载: 导出CSV
-
李鸿晶, 王通, 2011a.结构地震反应分析的逐步微分积分方法.力学学报, 43(2):430-435. http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=CAS201303040000629782 李鸿晶, 廖旭, 王通, 2011b.结构地震反应DQ解法的两种数值格式.应用基础与工程科学学报, 19(5):758-766. http://d.old.wanfangdata.com.cn/Periodical/yyjcygckxxb201105008 廖旭, 李鸿晶, 孙广俊, 2013.基于DQ原理的结构弹塑性地震反应分析.工程力学, 30(7):161-166. http://www.cnki.com.cn/Article/CJFDTotal-GCLX201307026.htm 王鑫伟, 1995.微分求积法在结构力学中的应用.力学进展, 25(2):232-240. doi: 10.3321/j.issn:1000-0992.1995.02.010 Bellman R., Casti J., 1971. Differential quadrature and long-term integration. Journal of Mathematical Analysis and Applications, 34 (2):235-238. doi: 10.1016/0022-247X(71)90110-7 Bellman R., Kashef B. G., Casti J., 1972. Differential quadrature:a technique for the rapid solution of nonlinear partial differential equations. Journal of Computational Physics, 10 (1):40-52. doi: 10.1016/0021-9991(72)90089-7 Bert C. W., Jang S. K., Striz A. G., 1988. Two new approximate methods for analyzing free vibration of structural components. AIAA Journal, 26 (5):612-618. doi: 10.2514/3.9941 Bert C. W., Wang X. W., Striz A. G., 1993. Differential quadrature for static and free vibration analyses of anisotropic plates. International Journal of Solids and Structures, 30 (13):1737-1744. doi: 10.1016/0020-7683(93)90230-5 Bert C. W., Malik M., 1996. The differential quadrature method for irregular domains and application to plate vibration. International Journal of Mechanical Sciences, 38 (6):589-606. doi: 10.1016/S0020-7403(96)80003-8 Bert C. W., Malik M., 1997. Differential quadrature:a powerful new technique for analysis of composite structures. Composite Structures, 39 (3-4):179-189. doi: 10.1016/S0263-8223(97)00112-8 Fung T. C., 2001a. Solving initial value problems by differential quadrature method——part 1:first-order equations. International Journal for Numerical Methods in Engineering, 50 (6):1411-1427. doi: 10.1002/(ISSN)1097-0207 Fung T. C., 2001b. Solving initial value problems by differential quadrature method——part 2:second-and higher-order equations. International Journal for Numerical Methods in Engineering, 50 (6):1429-1454. doi: 10.1002/(ISSN)1097-0207 Jang S. K., Bert C. W., Striz A. G., 1989. Application of differential quadrature to static analysis of structural components. International Journal for Numerical Methods in Engineering, 28 (3):561-577. doi: 10.1002/(ISSN)1097-0207 Kang K., Bert C. W., Striz A. G., 1995. Vibration analysis of shear deformable circular arches by the differential quadrature method. Journal of Sound and Vibration, 183 (2):353-360. http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=cd5d834cd6a527d56839a5876ae0d87a Kang K., Bert C. W., Striz A. G., 1996. Vibration analysis of horizontally curved beams with warping using DQM. Journal of Structural Engineering, 122 (6):657-662. doi: 10.1061/(ASCE)0733-9445(1996)122:6(657) Liew K. M., Han J. B., Xiao Z. M., 1996a. Differential quadrature method for thick symmetric cross-ply laminates with first-order shear flexibility. International Journal of Solids and Structures, 33 (18):2647-2658. doi: 10.1016/0020-7683(95)00174-3 Liew K. M., Han J. B., Xiao Z. M., et al., 1996b. Differential quadrature method for Mindlin plates on Winkler foundations. International Journal of Mechanical Sciences, 38 (4):405-421. doi: 10.1016/0020-7403(95)00062-3 Liu J., Wang X. W., 2008. An assessment of the differential quadrature time integration scheme for nonlinear dynamic equations. Journal of Sound and Vibration, 314 (1-2):246-253. doi: 10.1016/j.jsv.2008.01.004 Newmark N. M., 1959. A method of computation for structural dynamics. Journal of the Engineering Mechanics Division, 85 (3):67-94. http://cn.bing.com/academic/profile?id=04737dd237295d8ef1e115cb61eed609&encoded=0&v=paper_preview&mkt=zh-cn Sherbourne A. N., Pandey M. D., 1991. Differential quadrature method in the buckling analysis of beams and composite plates. Computers & Structures, 40 (4):903-913. doi: 10.1016-0045-7949(91)90320-L/ Striz A. G., Jang S. K., Bert C. W., 1988. Nonlinear bending analysis of thin circular plates by differential quadrature. Thin-Walled Structures, 6 (1):51-62. doi: 10.1016/0263-8231(88)90025-0 Wang X., Bert C. W., Striz A. G., 1993. Differential quadrature analysis of deflection, buckling, and free vibration of beams and rectangular plates. Computers & Structures, 48 (3):473-479. http://cn.bing.com/academic/profile?id=5a95130e65486d67576bb93101054491&encoded=0&v=paper_preview&mkt=zh-cn Wang X., Striz A. G., Bert C. W., 1994. Buckling and vibration analysis of skew plates by the differential quadrature method. AIAA Journal, 32 (4):886-889. doi: 10.2514/3.12071 Zeng H., Bert C. W., 2001. A differential quadrature analysis of vibration for rectangular stiffened plates. Journal of Sound and Vibration, 241 (2):247-252. doi: 10.1006/jsvi.2000.3295 期刊类型引用(0)
其他类型引用(1)
-
下载:
点击查看大图
计量
- 文章访问数: 96
- HTML全文浏览量: 96
- PDF下载量: 9
- 被引次数: 1
