《水工结构的三维阶谱有限元分析》由会员分享，可在线阅读，更多相关《水工结构的三维阶谱有限元分析（7页珍藏版）》请在金锄头文库上搜索。
1、水工结构的三维阶谱有限元分析程昭陈胜宏(武汉水利电力大学水电学院)摘要根据P型有限单元法的阶谱特点,详细论述了三维阶谱单元法的基本的分析过程和具体的 实现路径，包括基函数的构造、边界约束条件的处理、刚度矩阵和荷载列阵的形成、提高数值积 分效率的途径等。引入了三维阶谱单元法的虚结点和广义结点的概念，并把三维阶谱单元法应用 于水工结构计算。关键词有限元法，阶谱单元，三维分析，P型，水工结构。本文于1999年5月24日收到，系国家自然科学基金资助项0(59979021).有限元法是研究工程结构问题最为广泛的数值方法。从逼近真实解的途径分 类，有限兀法可以分为：(l)h型，即传统的有限元法，通过减小单
2、元尺寸h来 提高有限兀解的精度；(2)p型，通过增加基函数多项式的阶数p来提高有限元 解的精度;hp型，综合了(1)(2)两种方法。文献1给出了 p型有限元法的 理论分析，指出P型有限元的收敛速度比h型快，而对于奇点问题，至少是h 型的两倍。文献2阐述了 p型有限元法小的阶谱概念及其优点。基函数具有阶 谱特性的P型有限兀法，称为阶谱有限单兀法或阶谱单兀法。当前,P型有限元 法一般都采用阶谱单兀，所以在不作特别说明的情况下，P型有限兀法和阶谱单 元法是同一个概念。文献引对一维阶谱单元法作了较为详细的研究。FI前对阶 谱有限单元法的研究主要是针对二维问题，且偏重于理论分析，对工程计算的具 体过程如
3、何实现论述较少。比如，约束如何体现，如何提高数值积分效率等等。 本文详细地给出了三维阶谱单兀法的基本的分析过程和实现路径，包括基函数的 构造、边界约束条件的处理、刚度矩阵和荷载列阵的形成、提高数值积分效率的 途径等。引入了三维阶谱单元法中的虚结点和广义结点的概念，并应用于水工结 构工程计算，取得了满意的效果。1三维阶谱单元的基函数11基函数的特点 以比表示单元尺寸不变时pi阶阶谱单元逼近空间，阶谱的 要领即低阶单兀逼近空间是高阶单兀逼近空间的一个子集：H# UH旳UUH鸽U旧必U，O1去2 -Pi ） = 8+ 12max（0,p T） + 6max 0,- 2）（j _3） + max如果不
4、限制面基函数的i+j二p和体基函数的i+j+Ep,而分别改为i, jWp 和i, j, kWp即相当于Lagrange族基函数。Serendity族基函数比Lagrange族 基函数个数少，能减轻计算量，且性能很好，本文以前者为例作分析。其实，卩可以是一切在一维阶谱单元中可行的阶谱基函数。如果让上面的 p(t) = (l-t2) J(2, 2)p(t),其中 Jpi是 p 阶 Jacobian 单位(2, 2)多项式，就可以构成基于Jacobian多项式的Serendity族和Lagrange族阶谱基函数。2广义结点和边界约束条件三维阶谱单元法的拓扑信息包括：Q)各结点的坐标；(b)各个单元的
5、点编号 信息;(c)各个棱的点编号信息；(d)各个面的棱编号信息；(e)各个单元的面编 号信息。其中(c)-(e)可由(a)和(b)通过程序生成。由于阶谱单元法即使是粗糙 的网格也能达到很高的精度，不像传统的有限兀法需要精细地剖分网格，能大大 减小前处理的工作量。传统有限元法把面和棱的约束转化为面和棱上点的约束。由于阶谱单元法中 面和棱有各自的基函数，因此对单元贡献了自由度，可以看成是虚结点的自由度, 这些虚结点的自由度和实结点的自由度同等对待，虚结点和实结点可以统称为广 义结点。这样，点基函数、棱基函数、面基函数和体基函数就统一为广义结点的 基函数。有了广义结点的概念，许多问题就很容易理解。
6、例如，固定面的约束体 现为该面的四个角点为固定约束，该面的棱所代表的虚结点为固定约朿，该面所 代表的虚结点也为固定约束；其它约束情况可以此类推。因为基函数的阶数p 小于4时没有而基函数，p小于2时没有棱基函数，故当p小于4时不考虑而的 约束，当p小于2吋不考虑棱的约束。3刚度矩阵与荷载移植3.1单元刚度矩阵p二1吋，点18分别作为第18个结点，相应的基函数为 Nx-s； p=2时,增加棱112作为第920个结点(虚结点,下面省略),相应的 基函数为EiJ； p二3 11寸，增加棱112作为第2132个结点，相应的基函数为 Ej; P二4吋，增加棱112作为第3344个结点,相应的基函数为E_1
7、2,再 增加面16作为第4450个结点,相应的基函数为 计;P二5吋，增加棱1 12作为第51-62个结点，相应的基函数为-2,再增加面16作为第63 74个结点,相应的基函数为计一和计7卩二6时，增加棱112作为第75 86个结点，相应的基函数为-2,再增加而16作为第87104个结点，相 应的基函数为F巴十计3)_和计，再增加一个体结点作为第105个结点, 相应的基函数为22)；如此类推。基函数可以统一记为4,表示第i个广义 结点的基函数,其中lWiWfe(p).2-1网吗少=Z坷匕，W =2-1坐标的插值仍同等参单元法：8882-1ZH于是心=J耳D巧1 1 1 i打严D弓伸必必这样，单
8、元刚度矩阵就是阶谱形式。3.2整体自由度的整和设结构的单元数、实结点个数、棱的个数和面的个数分 别为N、汕、凡、N，则总的结点的个数为：fsp) = N氓 +max(O, - 1)A + maxM20,|(p-2Xp-3)眄+max 0,2(上_4；(1；将结构看作点数、棱数、而数和体数分别为陥汕，Nf, N的一个大的单元， 将汕，Ne, Nf, N分别对应3.1节中单元的点数8,棱数12,面数6和体数1,以 3.1节相同的方法编排整体结点，并去掉约束所对应的行和列，能使总体刚度矩 阵也成为阶谱形式。33荷载移植 荷载移植与传统的有限元法类似：r xI(11)k)=尸個氓=出心种矽恋(i =
9、1,2,/(刃)需要注意的是，由于虚结点的存在，移植后的结点荷载和与原荷载不存在平衡关 系，只是实结点的移植荷载和与原荷载平衡。这并不妨碍问题的解决。经3. 2 节方法对整体自由度进行整和后，整体荷载列阵F也是阶谱形式。4数值积分半阶数升高时，基函数将是高阶多项式，为保证精度，需耍更多的高斯积分 点，数值积分的运算量将随着阶数升高成倍增加，对三维问题尤为突出。利用阶 谱单元的特点，减少运算量的措施有：充分利用B的稀疏性。对于弹性问题，还可以利用D的稀疏性，将6 X6的矩阵B.J的各个元素先算出来，再对各个元素作数值积分。(2)单元刚度矩阵的不同元素采用不同阶高斯积分来计算。实际上，如果对 高阶
10、阶谱单元，各个兀素全部都用高阶高斯积分，计算吋间简直难以接受。G)另外从基函数的构造可以看到，即使一个方向上的阶数很高，另两个方 向上阶数可能很低，因此不同方向上也可以采用不同阶的高斯积分：dr ds4丿GMdr ds(12)其中,r, se 设P(f, )表示多项式函数f在方向上的阶数,II/2drdr(13)/2T才七4 6/宀八-7 表1单刚形成时间比较(单位:S)表1对不同积分方案的效率进行了P=1P=2p=3p二 4p=5方法(2)0. 4259. 73575. 495 644. 4955140. 450方法(3)0. 4254. 93024. 101 101.845 467. 12
11、8注：使用II 233微机，一般编程技巧。比较，从表1可以看出，采用措施能 在(2)的基础上进一步大大减少运算吋 间，当阶数较高时，甚至能减少90%以上 的运算时间。(4)采用 Howard E. Hinnant 提出的 矢量积分法叫也能大大减少计算时间, 这里不再赘述。5应用举例为便于分析整理成果，取段面如图1的混凝土重力坝，坝高100m,坝顶宽 10m,上游坡而垂直，下游坡而斜率1 ：0.75,坝体容重2400kg/m3,弹性模量 25GPa,泊松比1/6,坝基容重2600kg/m3,弹性模量75GPa,泊松比0.25.坝基 上游取1.5倍的坝高，下游取2倍的坝高，坝基的深度取2倍的坝高。
12、断面上的 单元网格如图1,再取10ni的坝段宽构成直四棱柱形式的三维六面体网格，-共 114个单元，284个实结点。以坝踵为坐标原点，x轴指向下游，z轴垂直向上， y轴由右手法则确定。两侧坝而无y向位移，坝基上下游而无x向位移，坝基底 面无z向位移。假设坝基地应力场由坝基重力场产生，坝体一次浇筑完成，水库 一次库满，下游无水。图1单元网格的侧视(实际是空间网格)控制总能量范数误差小于10%,实际计算当p=3时即达到要求。p=l, 2,3 的主要数据列表比较如表2,从表1、图2和图3可以看出，整体自由度随着阶 数的升高迅速增大，能量范数收敛较快，P二2和时位移和能量范数已相当接 近；坝踵处的应力等值线随着阶数的升高变得更密；不同的阶数情况下，应力的 计算结果差别最大的地方是在坝踵处，体现了坝踵处应力集中的特点。表2主要数据比较总刚密度最大位移/mm能量范数P二 182845105. 786%5. 992. 117E+1327. 04%P二 22093618754. 049%6. 742.135E+1316. 69%P二 332158832403. 795%6. 802.141E+139.93%(a)p=l时第一主应力等值线(b) p=2时第一主应力等值线
