在上篇几何非线性| 应变张量 ，得到拉格朗日应变表达式为
\mathbf E = \frac{1}{2}(\mathbf H + \mathbf H^T + \mathbf H^T \mathbf H)
用指标记法
E_{ij}=\frac{1}{2}(\frac{\partial u_i }{\partial X_j}+\frac{\partial u_j }{\partial X_i}+
\frac{\partial u_k }{\partial X_j}\frac{\partial u_k }{\partial X_i}) \quad (1)
对于杆系结构，有
\epsilon_{x}=\frac{\partial u }{\partial x}+\frac{1}{2}(\frac{\partial u }{\partial x})^2+
\frac{1}{2}(\frac{\partial v }{\partial x})^2+\frac{1}{2}(\frac{\partial w }{\partial x})^2 \quad (2)
拉格朗日应变适用于描述几何非线性。
▲图1
如图1所示的桁架单元，局部坐标下的位移插值
\begin{split}
u(x) &=[1- \frac{x}{l},0,\frac{x}{l},0]\begin{Bmatrix}
u_1 \\
v_1 \\
u_2 \\
v_2 \\
\end{Bmatrix} \\
v(x) &=[0,1- \frac{x}{l},0,\frac{x}{l}]\begin{Bmatrix}
u_1 \\
v_1 \\
u_2 \\
v_2 \\
\end{Bmatrix}
\end{split}
\mathbf u = \mathbf N\mathbf q^e \quad (3)
其中，
\mathbf q^e 是单元节点位移矩阵。
\begin{split}
\frac{\partial u }{\partial x} &=u^{'}\\
&= \frac{1}{l}\begin{bmatrix}
-1 & 0 & 1 & 0 \\
\end{bmatrix}
\begin{Bmatrix}
u_1 \\
v_1 \\
u_2 \\
v_2 \\
\end{Bmatrix} \\
&= \mathbf C \mathbf q^e
\end{split}
\begin{split}
\frac{\partial v }{\partial x} &=v^{'}\\
&= \frac{1}{l}\begin{bmatrix}
0 & -1 & 0 & 1 \\
\end{bmatrix}
\begin{Bmatrix}
u_1 \\
v_1 \\
u_2 \\
v_2 \\
\end{Bmatrix} \\
&= \mathbf D \mathbf q^e
\end{split}
拉格朗日应变
\epsilon = \mathbf C \mathbf q^e + \frac{1}{2} {\mathbf q^e}^T {\mathbf C}^T \mathbf C \mathbf q^e +\frac{1}{2} {\mathbf q^e}^T {\mathbf D}^T \mathbf D \mathbf q^e \quad (4)
虚位移
\delta \mathbf u = \mathbf N \delta \mathbf q^e \quad (5)
虚应变
\begin{split}
\delta \epsilon &= \frac{\partial \epsilon }{\partial \mathbf q^e }\delta \mathbf q^e \\
&= \mathbf B(\mathbf q^e) \delta \mathbf q^e
\end{split} \quad (6)
这里，(6)用到了变分运算公式
\delta f = \frac{\partial f }{\partial x }\delta x + \frac{\partial f }{\partial y }\delta y
内力虚功为
\begin{split}
\delta W_i & = \int_V \delta \boldsymbol {\epsilon}^T \boldsymbol{\sigma} dV \\
& = \delta \mathbf q^{eT} \int_V \mathbf B^T \boldsymbol
{\sigma} dV
\end{split} \quad (7)
记
f_i= \int_V \mathbf B^T \boldsymbol
{\sigma} dV \quad (8)
则
\delta W_i = \delta \mathbf q^{eT} f_i \quad (9)
由(6)可得
\begin{split}
\bf B &= \frac{\partial \epsilon }{\partial \mathbf q^{e}}\\
&= \frac{\partial \epsilon }{\partial u^{'}}\frac{\partial u^{'} }{\partial \mathbf q^{e}}+ \frac{\partial \epsilon }{\partial v^{'}}\frac{\partial v^{'} }{\partial \mathbf q^{e}}\\
&= (1+u^{'}) \mathbf {C} + v^{'} \mathbf D \\
&= \mathbf C + \mathbf q^{eT} \mathbf C^T\mathbf C + \mathbf q^{eT} \mathbf D^T\mathbf D
\end{split} \quad (10)
应力
\begin{split}
\sigma &= E(u^{'}+\frac{1}{2}u^{'2} + \frac{1}{2}v^{'2})\\
&= E(\mathbf C \mathbf q^e + \frac{1}{2} {\mathbf q^e}^T {\mathbf C}^T \mathbf C \mathbf q^e +\frac{1}{2} {\mathbf q^e}^T {\mathbf D}^T \mathbf D \mathbf q^e)
\end{split} \quad (11)
由(8)(10)(11)可得
f_i= ((1+u^{'}) \mathbf C^T + v^{'} \mathbf D^T )\sigma Al \quad (12)
▲图2
如图2所示的非线性迭代过程，当某一迭代步
i 达到收敛标准时，可以认为处于平衡状态，即
f_i(\mathbf q^i) = f_e(\mathbf q^i) \quad (13)
式中
f_i 是结构内力，
f_e 是外荷载，
\mathbf q^i 是
i 迭代步时的节点位移。
i+1 迭代步时的内力用一阶泰勒展开
f_i(\mathbf q^{i+1}) \approx f_i(\mathbf q^i) + \frac{\partial f_i }{\partial \mathbf q} \Delta \mathbf q^{i} \quad (14)
由(13)(14)得
\frac{\partial f_i }{\partial \mathbf q} \Delta \mathbf q^{i} = f_e(\mathbf q^{i+1})-f_i(\mathbf u^i) \quad (15)
记
\mathbf K_T = \frac{\partial f_i }{\partial \mathbf q} \quad (16)
其中，
\mathbf K_T 叫做切线刚度矩阵，(15)可写成
\mathbf K_T \Delta \mathbf q^{i} = f_e(\mathbf q^{i+1})-f_i(\mathbf q^i) \quad (17)
\mathbf K_T 是内力的导数，
f_e(\mathbf q^{i+1}) 是新的荷载步下的外荷载。
\begin{split}
\mathbf K_T &= \int_V \mathbf B^T \frac{\partial \boldsymbol
{\sigma} }{\partial \boldsymbol {\epsilon}}
\frac{\partial \boldsymbol {\epsilon}} {\partial \mathbf q}
dV + \int_V \frac{\partial \mathbf B^T }{\partial \mathbf q} \boldsymbol {\sigma}dV\\
&= \int_V \mathbf B^T E \mathbf B dV + \int_V \frac{\partial \mathbf B^T }{\partial \mathbf q} \boldsymbol {\sigma}dV\\
&= \mathbf K_{\mathbf q} + \mathbf K_{\boldsymbol
{\sigma}}\\
\end{split} \quad (17)
其中
\mathbf K_{\mathbf q} 叫做初始刚度矩阵，
\mathbf K_{\sigma} 叫做几何刚度矩阵。对于桁架单元
\begin{split}
\mathbf K_{\mathbf q} &= [(1+u^{'}) \mathbf {C}^T + v^{'} \mathbf D^T][(1+u^{'}) \mathbf {C} + v^{'} \mathbf D]EAl\\
&= [(1+u^{'})^2\mathbf {C}^T\mathbf {C}+ v^{'}(1+u^{'})(\mathbf {D}^T\mathbf {C}+\mathbf {C}^T\mathbf {D} ) + v^{'2}\mathbf {D}^T\mathbf {D}] EAl\\
\end{split} \quad (18)
几何刚度矩阵
\begin{split}
\mathbf K_{\boldsymbol
{\sigma}} &= \frac{\partial }{\partial u^{'} }(1+u^{'}) \mathbf {C}^T \frac{\partial u^{'} }{\partial \mathbf q } +\frac{\partial }{\partial v^{'} }(v^{'}) \mathbf {D}^T \frac{\partial v^{'} }{\partial \mathbf q } \boldsymbol
{\sigma}Al \\
&= (\mathbf {C}^T\mathbf {C}+\mathbf {D}^T\mathbf {D}) \boldsymbol
{\sigma}Al
\end{split} \quad (19)
其中
\begin{split}
\mathbf {C}^T\mathbf {C} &= \frac{1}{l^2}
\begin{bmatrix}
-1 \\
0 \\
1 \\
0 \\
\end{bmatrix}
\begin{bmatrix}
-1 & 0 & 1 & 0 \\
\end{bmatrix} \\
&= \frac{1}{l^2}
\begin{bmatrix}
1 & 0 & -1 & 0 \\
0 & 0 & 0 & 0 \\
-1 & 0 & 1 & 0 \\
0 & 0 & 0 & 0 \\
\end{bmatrix} \\
\end{split}
同理
\begin{split}
\mathbf {D}^T\mathbf {D} &= \frac{1}{l^2}
\begin{bmatrix}
0 & 0 & 0 & 0 \\
0 & 1 & 0 & -1 \\
0 & 0 & 0 & 0 \\
0 & -1 & 0 & 1 \\
\end{bmatrix} \\
\end{split}
最终得到局部坐标下的切线刚度矩阵为
\begin{split}
\mathbf K_T &= \frac{EA}{l}
\begin{bmatrix}
(1+u^{'})^2 & (v^{'}+v^{'}u^{'}) & -(1+u^{'})^2 & -(v^{'}+v^{'}u^{'}) \\
(v^{'}+v^{'}u^{'}) & v^{'2} & -(v^{'}+v^{'}u^{'}) & -v^{'2} \\
-(1+u^{'})^2 & -(v^{'}+v^{'}u^{'}) & (1+u^{'})^2 & (v^{'}+v^{'}u^{'}) \\
-(v^{'}+v^{'}u^{'}) & -v^{'2} & (v^{'}+v^{'}u^{'}) & v^{'2} \\
\end{bmatrix} \\
&\quad + \frac{\sigma A}{l}
\begin{bmatrix}
1 & 0 & -1 & 0 \\
0 & 1 & 0 & -1 \\
-1 & 0 & 1 & 0 \\
0 & -1 & 0 & 1 \\
\end{bmatrix} \\
\end{split} \quad (20)
局部坐标和整体坐标下的节点位移转换关系
\mathbf q^e = \mathbf T \mathbf q^g
整体坐标下的切线刚度矩阵
\mathbf K_T^g = \mathbf T^T \mathbf K_T\mathbf T
粤公网安备44030502008569号