p-细化(提高多项式阶数)
理论:加密还是提高阶数?
概述
网格收敛一般是指单元加密,但听说也可以提高单元阶数。哪个更好?
取决于解的光滑程度。加密单元的h-细化,误差只按单元尺寸的p次方减小。保持单元不变、提高形函数阶数的p-细化,解光滑时误差随自由度呈指数减小。但在尖角、集中载荷点等解不光滑处,提高阶数误差也难以减小。于是有两者结合的hp法。我们用一维算例在同样自由度下比较三者。
收敛速度
$\| \cdot \|_E$ 为能量范数(对应应变能的误差),$N$ 为自由度,$\lambda$ 为奇异性强度(解为 $x^{\alpha}$ 时 $\lambda = \alpha - 1/2$)。有奇异性时,h版和p版都受 $\lambda$ 限制。
层次形函数
p版常用积分Legendre多项式作为形函数。提高阶数时低阶函数保持不变,只追加新函数,刚度矩阵的计算可以复用。本算例使用这种形函数。
培育p版的人们
p版有限元法由美国圣路易斯华盛顿大学的Szabó等人于1970年代为结构分析开发,1981年Babuška、Szabó和Katz给出了收敛理论。光滑解指数收敛、有奇异点时通过设计网格也能指数收敛的hp法理论,由Babuška和Guo在1980年代给出。Szabó等人的研究也发展成保证应力集中部位精度的商用软件。
算例
例1:光滑解
一维问题 −u″ = f,u(0) = 0,解 u = sin(πx/2)。在自由度16下比较能量范数相对误差:
| 方法 | 相对误差 |
|---|---|
| h:16个一次单元 | 2.8×10⁻² |
| h:8个二次单元 | 1.4×10⁻³ |
| p:2个单元、阶数8 | 6.7×10⁻¹¹ |
保持2个单元把阶数从1提高到8,误差为 2.2×10⁻¹ → 2.3×10⁻² → 1.5×10⁻³ → 7.5×10⁻⁵ → 3.0×10⁻⁶ → 9.7×10⁻⁸ → 2.7×10⁻⁹ → 6.7×10⁻¹¹,阶数每提高1下降一个数量级以上。一次单元自由度加倍误差只减半。
例2:奇异解
同一问题,解取 u = x^0.7。x = 0处导数为无穷大,与尖角附近的应力性质相同。
| 方法 | 自由度 | 相对误差 |
|---|---|---|
| h:64个一次单元 | 64 | 1.9×10⁻¹ |
| h:64个二次单元 | 128 | 1.4×10⁻¹ |
| p:2个单元、阶数8 | 16 | 1.6×10⁻¹ |
| hp:比0.15几何分割7层 | 36 | 3.0×10⁻² |
有奇异点时,无论64个单元还是阶数8,误差都停在10%以上。
x^0.7 的 λ 只有0.2,均匀加密时单元尺寸每减半误差只减小13%。提高阶数也无法用光滑多项式表示尖锐的解。hp法向x = 0每层把单元尺寸缩小到0.15倍,与奇异点相接的单元保持一次,越往外阶数越高。每增加一层误差约变为0.68倍,自由度36时为3.0×10⁻²,约为均匀二次单元自由度128时的五分之一。实际结构分析中在应力集中部周围逐级布置小单元,也是同一思路。
选用要点
- 几何简单、解光滑的部分(薄壁弯曲、振动模态)p版有利。
- 在尖角、裂纹、接触端、材料界面附近逐级缩小单元。
- p版提高阶数后值仍持续变化时,怀疑有奇异点。
- 提高阶数时值的变化可以直接作为误差估计(p版收敛检查)。
- 高阶单元的几何(曲边界)也要用高阶表示,否则提高阶数没有效果。
“阶数提高应力也不收敛”
某设计者用p版分析软件检查L形板的应力,阶数从4提高到8,每次内角的最大应力都增加约10%。p版在光滑处几个阶数内值就稳定,不稳定本身就是“这里是奇异点”的信号。查图纸发现内角有半径2 mm的圆角。把圆角建进模型后,阶数6与8的差在0.5%以下,得到了可用于疲劳评估的值。
常见错误
错误与对策
| 错误 | 后果 | 对策 |
|---|---|---|
| 在奇异点只依赖p版 | 误差不减小 | 与几何分割结合 |
| 几何用一次近似 | 高阶效果消失 | 边界也用高阶表示 |
| 在畸变单元上高阶化 | 精度下降 | 改善单元形状 |
| 阶数过高 | 条件数恶化 | 使用层次形函数 |
| 只用局部应力判断 | 误判收敛 | 也检查能量 |
我想了解相关内容。