有限元动力学分析方程及解法VIP

本文格式为Word版，下载可任意编辑——有限元动力学分析方程及解法动力分析中平衡方程组的解法

1前言

描述结构动力学特征的基本力学变量和方程与静力问题类似，但所有的变量都是时间的函数。

基本变量

三大类变量ui(?,t)、?ij(?,t)和?ij(?,t)是坐标位置?(x,y,z)和时间t的函数，一般将其记为ui(t)?ij(t)?ij(t)。

基本方程(1)平衡方程

利用达朗贝尔原理将惯性力和阻尼力等效到静力平衡方程中，有

??i(t)??u?i(t)?0(1)?ij,j(t)?bi(t)??u

其中?为密度，?为阻尼系数。(2)几何方程

12?ij(t)?(ui,j(t)?uj,i(t))(2)

(3)物理方程

?ij(t)?Dijkl?kl(t)(3)

其中Dijkl为弹性系数矩阵。(4)边界条件

位移边界条件BC(u)为，

ui(t)?ui(t)在Su上(4)

力的边界条件BC(p)为，

?ij(t)nj?pi(t)在Sp上(5)

初始条件

ui(?,t?0)?ui0(?)(6)?i0(?)(7)?i(?,t?0)?uu虚功原理

基于上述基本方程，可以写出平衡方程及力边界条件下的等效积分形式，

??i??u?i?b)?uid???(?ijnj?pi)dA?0(8)?????(?ij,j??u?Sp对该方程右端第一项进行分部积分，并应用高斯-格林公式，整理得，

???(Dijklij??i?ui??u?i?ui)d??(?bi?uid???pi?uidA)?0(9)???kl??u?Sp有限元分析列式单元的节点位移列阵为，

Ute(t)?[u1(t),v1(t),w1(t),u2(t),v2(t),w2(t)?uk(t),vk(t),wk(t)](10)

单元内的插值函数为，

u(?,t)?N(?)Ute(t)(11)

其中N(?)为单元的形状函数矩阵，与相应的静力问题单元的形状函数矩阵完全一致，?为单元中的几何位置坐标。

基于上面的几何方程和物理方程及(11)式，将相关的物理量表达为节点位移的关系，有，

?(?,t)?[?]u(?,t)?[?]N(?)Ute(t)?B(?)Ute(t)(12)?(?,t)?D??DB(?)Ute(t)?S(?)Ute(t)(13)

?e(t)(14)?(?,t)?N(?)Uut??e(t)(15)??(?,t)?N(?)Uut将(12)-(15)供稿到虚功方程(9)中，有，

??e(t)?CeU?e(t)?KeUe(t)?Re(t)]T????[MeU?Ute(t)?0(16)tttt由于?Ute(t)具有任意性，消去该项并简写有，

??e?CeU?e?KUe?Re(17)Utttt其中，

Me??e??NTNd?(18)

Ce???NTNd?(19)

?eKe??e?BTDBd?(20)

Me为单元质量矩阵，Ce为单元阻尼矩阵，Ke为单元刚度矩阵。同样，将单元

的各个矩阵进行组装，可形成系统的整体有限元方程，即，

???CU??KU?R(21)MU??,U?其中M，C和K分别是系统的质量、阻尼和刚度矩阵，R是外荷载向量，U和U分别是有限元分割体的加速度、速度和位移向量。方程(21)是通过考虑在时刻t的静力平衡而推导出来的。

对静力或动力分析的选择(即在分析中是考虑或忽略与速度及加速度有关的力)，一般取决于工程上的判断，其目的在于减少所需要的分析工作量。但是，应当认识到，一个静力分析的假定，应当有理由说明它是正确的，否则，分析的结果就是无意义的。确实，在非线性分析中，采用忽略惯性力和阻尼力的假定，可能严重到难以求得甚至无法求得解答。

在数学上，方程(21)是一个二阶线性微分方程组，原则上可用求解常系数微分方程组的标准过程来求得方程组的解。但是，假使矩阵的阶数很高，则采用求解一般微分方程组的过程可能要付出很高的费用，除非特别利用系数矩阵K，C和M的特别性质。因此，在实用的有限元分析中，主要对几种有效的方法感兴趣，下面将集中介绍这几种方法。我们所考虑的基本过程，可分为两种求解方法：直接积分法和振型叠加法。初看起来，这两种方法似乎完全不同，但事实上它们有着密切的关系，至于选择这种或那种方法，只取决于它们的数值效果。

2直接积分法

在直接积分中对方程(21)是逐步地进行数值积分的，“直接〞的意思是，进行数值积分前没有进行把方程变为另一种形式的变换。实质上，直接积分是基于下面的两个想法，第一个想法是只在相隔?t的一些离散的时间区间上而不是试图在任一时刻t上满足方程(21)即包含有惯性力和阻尼力作用的(静力)平衡是在求解区间上的一些离散时刻点上获得的。因此，似乎在静力分析中使用过的所有求解方法，在直接积分法中或许也能有效地使用；其次个想法是假定位移、速度和加速度在每一时间区间?t内变化。

?,U??来表示初始时刻(t?0)的位移、速度和加速度向下面假设分别用U0,U00量为已知，要求出方程(21)从t?0到t?T的解。在求解时，把时间全程T划分为几个相等的时间区间?t(即?t?T/n)，所用的积分格式是在时刻0,?t，

2?t,?,t,t??t,?,T上确定方程的近似解。由于计算下一个时刻的解的算法要考

虑到前面各个时刻的解，因此假定在时刻0,?t,?,t的解为已知，来推导出求时刻t??t的解的算法。计算时刻t??t的解对于计算自此以后?t的时刻上的解是有代表意义的，这样就可建立用来计算在所有离散时间点上解的一般算法。

（a）中心差分法

若把式(21)的平衡关系看作是一个常系数常微分方程组，便可以用任一有限差分表达式通过位移来近似表示加速度和速度。因此，在理论上，大量不同的有限差分表达式均可使用。但是，我们要求求解格式必需是有效的，这样便只需考虑少数几种计算格式。对某些问题求解是十分有效的一个过程是中心差分法，这个方法假定

???1UUtt??t?2Ut?Ut??t2?t(22)

??1?UUtt??t?Ut??t2?t????将式(22)代入t时刻的式(21)，可得

121?1????1?C?Ut??t?Rt??K?2M?Ut??2M?C?Ut??t(23)?2M?2?t?2?t??t??t????t

从式(23)我们可以求出Ut??t。应当注意，Ut??t的解是基于利用在时刻t的平衡条件。因此，该积分过程称为显式积分方法，且这样的积分格式在逐步解法中不需要对(有效)刚度矩阵进行分解。另一方面，以后所考虑的Houbolt，Wilson?及Newmark方法，要利用在t??t上的平衡条件，因而称为隐式积分方法。

另外还应注意到，应用中心差分法时，Ut??t的计算包含有Ut和Ut??t，因此，

?,U??都是已知的，计算在时刻?t的解，必需用一个具体的起始过程。由于U0,U00由关系式(22)可求

U??t?t2????U0??tU0?U0(24)

2具体计算步骤为

A．初始计算1．2．

形成刚度矩阵K、质量矩阵M和阻尼矩阵C。

?,U??。计算初始值U0,U003．4．5．6．7．

选取时间步长?t，要求?t??tcr（临界值）。计算系数a0?111a?，，，。a?a?2a13202?t?t2a2??aU??计算U??t?U0??tU030。

??aM?aC。形成有效质量矩阵M01T??LDL?作三角分解：M对M

B．每一时间步长内的计算1．

计算在时刻t的有效荷载：

??R??K?aM?U??aM?aC?U。Rtt2t01t??t2．3．化为

1?(25)MUt??t?Rt2?tT?。计算时刻t??t的位移：LDLUt??t?Rt必要时，依照式（11.3）计算时刻t速度和加速度。

假设所考虑的系统没有物理阻尼，即C是零矩阵，在这种情形下式(23)可简

其中

??R??K?aM?U??aM?aC?URtt2t01t??t因此，假使质量矩阵是对角形的，则解方程组（11.1）时就不需要进行矩阵的分

?，从而利用解，即只需进行矩阵相乘便可求得右端项的有效荷载向量RtU(i)t??t??R(i)t??t2???m??(26)?ii?)?(i)分别