线性静态有限元分析是有限元分析的最基本内容,通过求解结构的平衡方程,可以求得结构的位移,进而得到结构的应变、应力与约束反力。第2、3、5、6章通过不同方法得到了由单元平衡方程组装而成的结构平衡方程,在静态有限元分析之前需要先求解积分形式的单元刚度矩阵与载荷列阵。通常被积函数的表达式项数很多,进行精确的解析积分十分困难,只能用数值积分方法近似求积。
7.1 数值积分方法
数值积分的基本思想是把积分化为求和,如下式所示:
\[\int_{\xi_1}^{\xi_2} F(\xi) d\xi = \sum_{i} \alpha_i F(\xi_i) \]\[\int_{\eta_1}^{\eta_2} \int_{\xi_1}^{\xi_2} F(\xi, \eta) d\xi d\eta = \sum_{i} \sum_{j} \alpha_{ij} F(\xi_i, \eta_j) \]一维问题标准形式可写成
\[\int_{-1}^{+1} F(\xi) d\xi \]注意积分上、下限分别为 +1,-1。
2. 高斯积分法
如果事先不确定取样点的位置,对每个取样点有两个待定的量,即积分点的位置\(\xi_i\)和权系数\(H_i\)。如果有\(n\)个取样点就有\(2n\)个条件,可以确定一个\(2n - 1\)阶多项式,用这个多项式近似被积函数并精确积分,积分可写成
\[\int_{-1}^{+1} F(\xi) d\xi=\sum_{i = 1}^{n} H_i F(\xi_i) \](7.8)
表7 - 1列出了前三阶高斯积分的积分点坐标和对应的权系数。显然,选取相同数目取样点,高斯给出更高的积分精度。在有限元法中几乎毫无例外地都用高斯法进行数值积分。
表7 - 1 高斯积分的积分点坐标和权系数
积分点数\(n\) | 积分点坐标\(\xi_i\) | 积分权系数\(H_i\) |
---|---|---|
1 | \(\xi_1 = 0.00000000\) | \(H_1 = 2.00000000\) |
2 | \(\xi_1, \xi_2=\pm 0.57735027\) | \(H_1, H_2 = 1.00000000\) |
3 | \(\xi_1, \xi_3=\pm 0.77459667\) | \(H_1, H_3 = 0.55555556\) |
\(\xi_2 = 0.00000000\) | \(H_2 = 0.88888889\) |
正因为构造的被积函数是\(2n - 1\)次多项式,因此\(n\)个积分点的高斯积分可达\(2n - 1\)阶的精度。也就是说如果\(F(\xi)\)是不超过\(2n - 1\)次的多项式,积分结果将是精确的。
7.2 单元刚度矩阵与载荷列阵的解析积分和数值积分
对于两节点杆单元,\(K^e\)可以解析积分出具体数值,而不必采用数值积分。将
将(4.16)式代入(6.21)式,得
\[K^e = EA \int_{x_1^e}^{x_2^e} \left( \frac{1}{l^e} \left[ \begin{array}{cc} -1 & 1 \\ 1 & -1 \end{array} \right] \right)^T \left( \frac{1}{l^e} \left[ \begin{array}{cc} -1 & 1 \\ 1 & -1 \end{array} \right] \right) dx \]\[= \frac{EA}{(l^e)^2} \int_{x_1^e}^{x_2^e} \left[ \begin{array}{cc} -1 & 1 \\ 1 & -1 \end{array} \right] \left[ \begin{array}{cc} -1 & 1 \\ 1 & -1 \end{array} \right] dx = \frac{EA}{l^e} \left[ \begin{array}{cc} 1 & -1 \\ -1 & 1 \end{array} \right] \]对于3节点2次杆单元(见4.3节),单元应变矩阵\(B^e\)为
\[B^e = \frac{dN^e(x)}{dx} = \left[ \begin{array}{ccc} \frac{dN_1^e(x)}{dx} & \frac{dN_2^e(x)}{dx} & \frac{dN_3^e(x)}{dx} \end{array} \right] = \frac{1}{(l^e)^2} \left[ \begin{array}{ccc} 2x - (x_2^e + x_3^e) & -4x + 2(x_2^e + x_3^e) & 2x - (x_2^e + x_3^e) \end{array} \right] \]令
\[x_c = \frac{x_2^e + x_3^e}{2} \]由于\(x_2^e - x_1^e = l^e\),\(x_3^e - x_1^e = x_3^e - x_2^e = l^e/2\),那么(7.10)式可进一步表达为
\[B^e = \frac{1}{(l^e)^2} \left[ \begin{array}{ccc} (x - x_c) - \frac{l^e}{2} & -2(x - x_c) & (x - x_c) + \frac{l^e}{2} \end{array} \right] \]将(7.12)式代入下式的单元刚度矩阵:
\[K^e = \int_{x_1^e}^{x_3^e} (B^e)^T A E B^e dx \]其中
\[(B^e)^T B^e = \left[ \begin{array}{ccc} \frac{4(x - x_c)^2}{(l^e)^4} & \frac{8(x - x_c)}{(l^e)^3} & \frac{32(x - x_c)^2}{(l^e)^4} \\ \frac{8(x - x_c)}{(l^e)^3} & \frac{32(x - x_c)^2}{(l^e)^4} & \frac{64(x - x_c)^3}{(l^e)^5} \\ \frac{16(x - x_c)^2}{(l^e)^4} & \frac{32(x - x_c)^3}{(l^e)^5} & \frac{1}{(l^e)^4} \end{array} \right] \]\[+ \left[ \begin{array}{ccc} \frac{16(x - x_c)^2}{(l^e)^4} & \frac{8(x - x_c)}{(l^e)^3} & \frac{32(x - x_c)^2}{(l^e)^4} \\ \frac{8(x - x_c)}{(l^e)^3} & \frac{32(x - x_c)^2}{(l^e)^4} & \frac{64(x - x_c)^3}{(l^e)^5} \\ \frac{4(x - x_c)^2}{(l^e)^4} & \frac{8(x - x_c)}{(l^e)^3} & \frac{1}{(l^e)^2} \end{array} \right] \]将(7.14)式代入(7.13)式并积分,得
\[K^e = AE \left[ \begin{array}{ccc} \frac{(l^e)^2}{12} \left( \frac{4(x - x_c)}{(l^e)^2} - \frac{1}{l^e} \right) & \frac{4x(x - 2x_c)}{(l^e)^3} - \frac{32(x - x_c)^3}{3(l^e)^4} \\ \frac{4x(x - 2x_c)}{(l^e)^3} - \frac{32(x - x_c)^3}{3(l^e)^4} & \frac{32(x - x_c)^3}{(l^e)^4} - \frac{64(x - x_c)^3}{(l^e)^5} \\ \frac{16(x - x_c)^3}{(l^e)^4} - \frac{4x(x - 2x_c)}{(l^e)^3} & \frac{32(x - x_c)^3}{(l^e)^4} - \frac{32(x - x_c)^3}{(l^e)^4} \\ \frac{16(x - x_c)^3}{(l^e)^4} - \frac{4x(x - 2x_c)}{(l^e)^3} & \frac{32(x - x_c)^3}{(l^e)^4} - \frac{32(x - x_c)^3}{(l^e)^4} \\ \frac{(l^e)^2}{12} \left( \frac{4(x - x_c)}{(l^e)^2} + \frac{1}{l^e} \right) \end{array} \right]_{x_1^e}^{x_3^e} \](7.15)式可进一步整理为
\[K^e = \frac{AE}{3l^e} \left[ \begin{array}{ccc} 7 & -8 & 1 \\ -8 & 16 & -8 \\ 1 & -8 & 7 \end{array} \right] \]显然,以上的解析积分是比较繁琐的,但是仍然可以手动推导,只是工作量较大。事实上,绝大多数的有限元单元都必须通过数值积分来计算。(7.14)式仅仅是3个自由度的杆单元,尚且如此复杂,而平面三角形单元、四边形单元等,其刚度矩阵\(K^e\)只能通过数值积分得到,而推导显式表达式的计算量将不可估计。
接下来,采用高斯数值积分求解(7.14)式,令
\[f(x) = EA \left[ \begin{array}{ccc} \frac{4}{(l^e)^2} (x - x_c) - \frac{1}{l^e} & \frac{8}{(l^e)^3} (x - x_c) & \frac{4}{(l^e)^2} (x - x_c) + \frac{1}{l^e} \\ \frac{8}{(l^e)^3} (x - x_c) & \frac{32}{(l^e)^4} (x - x_c)^2 & \frac{8}{(l^e)^3} (x - x_c) \\ \frac{4}{(l^e)^2} (x - x_c) + \frac{1}{l^e} & \frac{8}{(l^e)^3} (x - x_c) & \frac{4}{(l^e)^2} (x - x_c) + \frac{1}{l^e} \end{array} \right] \]为了将积分区间\([x_1^e, x_3^e]\)换到标准积分区间\([-1,1]\),可经变换
\[x = \frac{x_3^e - x_1^e}{2} \xi + \frac{x_3^e + x_1^e}{2} = \frac{l^e}{2} \xi + x_c, \quad x \in [x_1^e, x_3^e] \]且
\[dx = \frac{l^e}{2} d\xi \]此时
\[K^e = \int_{x_1^e}^{x_3^e} f(x) dx = \int_{-1}^{1} f\left( \frac{l^e}{2} \xi + x_c \right) \frac{l^e}{2} d\xi = \frac{l^e}{2} \int_{-1}^{1} g(\xi) d\xi \]其中,\(g(\xi) = f\left( \frac{l^e}{2} \xi + x_c \right) \frac{l^e}{2}\)。
首先采用单点高斯积分
\[K^e = \frac{l^e}{2} \int_{-1}^{1} g(\xi) d\xi = \frac{l^e}{2} H_1 g(\xi_1) \]其中
\[\xi_1 = 0, \quad H_1 = 2.0 \](7.22)
当\(\xi_1 = 0\)时
\[g(0) = \frac{EA}{(l^e)^2} \left[ \begin{array}{ccc} 1 & 0 & -1 \\ 0 & 0 & 0 \\ -1 & 0 & 1 \end{array} \right] \]将(7.22)式与(7.23)式代入(7.21)式,得到单点积分的刚度矩阵结果为
\[K^e = \frac{EA}{l^e} \left[ \begin{array}{ccc} 1 & 0 & -1 \\ 0 & 0 & 0 \\ -1 & 0 & 1 \end{array} \right] \]与精确解(7.16)相比较,可知采用单点积分误差较大。
(7.20)式的2点高斯积分公式为
\[K^e = \frac{l^e}{2} \int_{-1}^{1} g(\xi) d\xi = \frac{l^e}{2} [H_1 g(\xi_1) + H_2 g(\xi_2)] \]其中
\[\xi_1 = -0.57735027, \quad H_1 = 1 \]\[\xi_2 = +0.57735027, \quad H_2 = 1 \]对于\(\xi_1=-0.57735027\)
\[g(\xi_1)=\frac{EA}{(l^e)^2} \left[ \begin{array}{ccc} 4.6427344 & -4.9760678 & 0.3333333 \\ -4.9760678 & 5.3333333 & -0.3572655 \\ 0.3333333 & -0.3572655 & 0.0239226 \end{array} \right] \]对于\(\xi_2=+0.57735027\)
\[g(\xi_2)=\frac{EA}{(l^e)^2} \left[ \begin{array}{ccc} 0.0239226 & -0.3572655 & 0.3333333 \\ -0.3572655 & 5.3333333 & -4.9760678 \\ 0.3333333 & -4.9760678 & 4.6427344 \end{array} \right] \]将(7.27)式、(7.28)式代入(7.25)式,得到2点高斯积分的刚度矩阵为
\[K^e=\frac{EA}{l^e} \left[ \begin{array}{ccc} 2.3333333 & -2.6666667 & 0.3333333 \\ -2.6666667 & 5.3333333 & -2.6666667 \\ 0.3333333 & -2.6666667 & 2.3333333 \end{array} \right] \]由此可见,采用2点高斯积分已经可以得到精确结果。
(7.20)式的3点高斯积分公式为
\[K^e=\frac{l^e}{2} \int_{-1}^{1} g(\xi) d\xi=\frac{l^e}{2} [H_1 g(\xi_1)+H_2 g(\xi_2)+H_3 g(\xi_3)] \]其中
\[\xi_1=-0.77459667, \quad H_1=0.55555556 \]\[\xi_2=0, \quad H_2=0.88888889 \]\[\xi_3=+0.77459667, \quad H_3=0.55555556 \]对于\(\xi_1=-0.77459667\)
\[g(\xi_1)=\frac{EA}{(l^e)^2} \left[ \begin{array}{ccc} 6.49838658 & -7.89838653 & 1.39999994 \\ -7.89838653 & 9.59999977 & -1.70161325 \\ 1.39999994 & -1.70161325 & 0.30161330 \end{array} \right] \]对于\(\xi_2=0\)
\[g(\xi_2)=\frac{EA}{(l^e)^2} \left[ \begin{array}{ccc} 1 & 0 & -1 \\ 0 & 0 & 0 \\ -1 & 0 & 1 \end{array} \right] \]对于\(\xi_3 = +0.77459667\),
\[g(\xi_3)=\frac{EA}{(l^e)^2} \left[ \begin{array}{ccc} 0.30161330 & -1.70161341 & 1.40000007 \\ -1.70161341 & 9.60000027 & -7.89838653 \\ 1.40000007 & -7.89838653 & 6.49838658 \end{array} \right] \]将(7.32)式~(7.34)式代入(7.30)式,得到
\[K^e=\frac{EA}{l^e} \left[ \begin{array}{ccc} 2.33333333 & -2.66666667 & 0.33333333 \\ -2.66666667 & 5.33333333 & -2.66666667 \\ 0.33333333 & -2.66666667 & 2.33333333 \end{array} \right] \]便得到\(K^e\)的数值表达式,其与2点高斯积分结果一致。
由7.1节可知:当积分点的个数为\(n\)时,高斯积分具有\(2n - 1\)次的代数精度,即对于一切不高于\(2n - 1\)次的多项式,高斯积分都精确成立;若多项式次数高于\(2n - 1\)次,则不成立。那么,当积分点个数为1时,高斯积分具有1次代数精度;当积分点个数为2时,高斯积分具有3次代数精度。所以,对于3节点杆单元的2次多项式被积函数(7.17),需要2个积分点即可精确积分,这个结论也得到了上文的验证。
对于有3个及3个以上节点的单元,内部节点自由度可以在单元层次凝聚掉,而只保持外部节点自由度参加系统方程的集成,以提高计算效率。
之所以采用数值积分,是因为对于自由度较多的单元(如板壳结构单元、平面或实体单元)解析推导积分的过程太过繁琐,而并不是因为被积函数不可解析积分。事实上,有限元单元常用的Lagrange插值函数或Hermite插值函数都是解析可积的。
以上是单元刚度矩阵的积分,对于(6.22)式载荷列阵的积分相对容易,可以采用解析积分或者数值积分。此处,将2节点杆单元插值形函数代入(6.22)式的\(f_{x}^{e}\)项,同时假设\(b\)为常数,得到
\[f_{x}^{e}=\int_{x_1^{e}}^{x_2^{e}} (N^e)^T b dx = \frac{b l^e}{2} \left[ \begin{array}{c} 1 \\ 1 \end{array} \right] \]可以看出,(7.36)式的积分相当于将体力集中于杆单元的两个节点自由度上。另外,对于(6.22)式的\(f_{r}^{e}\)项,可由插值函数的基本属性\(N^e(x_r^{e})=\delta_{ij}\),得
\[f_{r}^{e}=(N^e)^T A \left. \frac{du^e}{dx} \right|_{x = x_r^{e}} = \left[ \begin{array}{c} 0 \\ A \left. \frac{du^e}{dx} \right|_{x = x_r^{e}} \end{array} \right] \]上式说明,作用于单元节点上的载荷,直接将其放置到单元载荷列阵相应自由度上即可。
标签:xi,frac,高斯,积分,right,array,left From: https://www.cnblogs.com/redufa/p/18641841