XII 线性屈曲分析
12.1 引言
本章主要由Markus Kriesch和André Wehr (Universitǎt der Bundeswehr München )编写，Matthias Goelke同时参与。
从2000年开始，由于商品价格的不断上涨，汽车对低油耗、低排放的要求越来越高。然而，追求轻量化设计（在一定程度上）抵消了由被动与主动安全系统和稳定增加的乘客舒适水平导致的附加重量。这就在许多行业内都带来了对智能轻量化概念的需求，而不仅仅是航空与汽车行业。
除了轻量化设计，正在增加使用优化软件分析对薄壁细长件进行轴向加载而引起的屈曲。
12.2 屈曲
细薄结构受到压缩载荷,还未达到材料强度极限而出现的失效状态称为屈曲。屈曲的特点是结构件在受到高压应力时突然失效，而失效点的实际压应力小于材料所能承受的极限压应力。
(http://en.wikipedia.org/wiki/Buckling)
换句话说，细长杆件（比如梁，其长度远大于其横截面）受到临界载荷时，会倒向一侧而不是继续承受更多载荷。
图: 细长梁力与位移的关系
(摘自Ihlenburg: Skript Technische Mechanik 2.7 Knicken. HAW Hamburg – FB MP; 2012)
此失效行为可用众所周知的线性屈曲分析来处理，目的是确定屈曲载荷因子λ和临界屈曲载荷。
在RADIOSS中，如果载荷因子λ > 1，则部件视为安全（比如，在屈曲产生前用实际载荷乘以λ）。
12.3 弹性屈曲
图: 欧拉屈曲案例; K=有效屈曲长度因子; L=柱的长度
(来自Läpele, Volker: Einführung in die Festigkeitslehre; Vieweg+Teubner Verlag; 2011)
1757年，Leonhard Euler 推导出以下方程：
Fcrit = π2 Ε Ι / (K L)2 = π2 Ε Ι / s2
在此
F = 最大或临界力（柱向的垂直力）
E = 弹性模量
I = 面积惯性矩（第二面积矩）
L = 柱长
K = 柱有效长度因子，其值取决于柱端的支撑状态，如下：
两端铰接（自由旋转），K = 1.0
两端固定，K = 0.50
一端固定，一端铰接，K = 0.707
一端固定，一端自由，K = 2.0
KL = 柱的有效屈曲长度
换句话说，临界力取决于：
• 柱长
• 横截面（第二表面积矩）
• 材料属性（杨氏模量，若为弹性材料）
• 边界条件（边界条件决定弯曲模态和变形的圆柱上各拐点的距离，拐点的距离越近，圆柱导致的结果越高） 圆柱的强度因材料分配而增加以提高惯性矩。可通过将材料尽可能分布到远离中心轴的横截面而不增加圆柱的重量来实现，同时保持材料厚度而不发生局部屈曲。这就验证了一个众所周知的事实，对于柱体来讲，一个管截面会比一个实心结构效果更好。
注：真实的建筑常常存在一些缺陷，比如，预变形，载荷尚未达到理想临界值而发生大变形或失效。线性屈曲分析常常高估了结构的强度与稳定性，从而导致非保守结果。因此，它不应作为唯一的衡量标准。然而，线性屈曲分析至少提供了预期的变形形状信息。
看看欧拉方程，并除以面积A定义屈曲发生的临界应力：
σkrit = Fcrit / A ≤ Re
其中Re为弹性极限：
σkrit = Fcrit / A = π2 E I / A s2 = π2 E / λ2 ≤ Re
且
λ = s / √ (I / A)
12.4 RADIOSS线性稳定性分析
在有限元分析中的线性屈曲分析，首先通过对结构施加一个参考水平的载荷线性Fref来解决。施加的是一个理想的单位载荷F。单位载荷和相应的约束，SPC，在第一个载荷步、工况下引用。 接下来执行一个标准的线性静态分析来获取用以生成几何刚度矩阵KG的应力。接着，其屈曲载荷会做为第二个载荷步、工况的一部分通过求解特征值问题计算出来。
(K-λKG)x=0
K是结构的刚度矩阵，λ是参考载荷的乘数。特征值问题的求解通常成生n个λ（屈曲载荷因子），其中n是自由度的数目（实际上，仅会计算出特征值的一个子集）。向量x是与特征值对应的特征向量。
特征值问题通过一个称之为“Lanczos法”的矩阵法来求解。并非所有的特征值都是必需的，对于屈曲分析只需要计算最小特征值的一小部分。
最小特征值与屈曲有关，临界或屈曲载荷为：
Fcrit = λcrit Fref
换句话说,
λcrit = Fcrit / Fref
因此
λc
λc > 1 safe
注：屈曲分析中所得到的位移结果可将屈曲模态的振型描述出来。
任何位移结果都是无意义的，同样适用于屈曲分析中的应力应变结果。
从理论到实践: 怎么样去设置屈曲分析
为了运行线性屈曲分析，需要以下两步：
第一步：
三个载荷集与两个载荷步、工况必须定义：
• 约束的载荷集 (SPC =单点约束; 不需要Card Image)
• 载荷的载荷集 (理想的单位载荷)；不需要Card Image
• 另一个载荷集(需要Card Image EIGRL)定义屈曲模态的数量(card image 如下)
• ND = 1: 使RADIOSS提取出第一阶屈曲模态
第二步
定义两个工况：一个静态和一个屈曲。载荷集SPCs和单位载荷定义一个静态工况。注意此工况的type设为linear static：
接下来定义屈曲工况（需要静态工况的结果）：
如前，SPC引用模型约束，STATSUB (BUCKLING)引用之前的静态工况，最后，METHOD引用提取屈曲模
态根数目信息的载荷集（例如EIGRL）。同时要注意，此时type要设为linear buckling。
注：
• STATSUB不可调用带有惯性释放的工况。
• 屈曲分析将忽略0维单元、MPC、RBE3和CBUSH单元等。
• 这些单元可用于屈曲分析，但对几何刚度矩阵KG无作用。
• 默认情况下，不包括刚性单元对几何刚度矩阵的贡献。用户必须在bulk data片段中添加PARAM、KGRGD、
YES使之包括刚性单元对几何刚度矩阵的贡献。
• 另外，用户可通过EXCLUDE工况信息条目，决定忽略其他单元对几何刚度矩阵的贡献，使用户有效的去
控制结构的哪些部分可用于屈曲分析。这此去除掉的属性仅仅是从几何刚度矩阵中移除，从而使屈曲分析处于弹性边界状况。这就意味着这些去除的属性依然会在屈曲模态中显示其运动。
例题：梁结构的线性屈曲分析
第一步
打开HyperMesh，选择RADIOSS (Bulk Data) 用户界面并读取文件BAR_EULER_1. Hm：
图: 导入模型“BAR_Euler_1.hm”. 1D单元显示为3D单元模式,
模型信息：
b = 10 mm
h = 20 mm
l = 1000 mm
E = 210.000 N/mm2 , ν=0.3 (steel)
面积惯性矩 = 1.666,667 mm4
分析屈曲临界力 Fcrit = 863,59 N (欧拉案例1)
模型用1D单元(类型为CBAR) 创建，其属性集为PBARL (若想了解1D单元的更多信息，请参考此书的相关章节). 材料参数为钢材。
为了更好地了解模型的详细信息，请回顾属性与材料各自的信息。
第二步
创建3个载荷集：
1. SPC (单点约束；无Card Image)
2. 力 静态单位力；无 Card Image)
3. 屈曲(使用卡片EIGRL)
载荷集很容易地通过右键单击Model Browser窗口，并激活选项Create > Load Collector创建。
接着，赋予各自的名称与颜色。
图：创建SPC和力 (不需
Card image)
图: 导入模型“BAR_Euler_1.hm”. 1D单元显示为3D单元模式,
模型信息：
b = 10 mm h = 20 mm l = 1000 mm
E = 210.000 N/mm2 , ν=0.3 (steel) 面积惯性矩 = 1.666,667 mm4
分析屈曲临界力 Fcrit = 863,59 N (欧拉案例1)
模型用1D单元(类型为CBAR) 创建，其属性集为PBARL (若想了解1D单元的更多信息，请参考此书的相关章节). 材料参数为钢材。
为了更好地了解模型的详细信息，请回顾属性与材料各自的信息。
第二步
创建3个载荷集：
1. SPC (单点约束；无Card Image) 2. 力 静态单位力；无 Card Image) 3. 屈曲(使用卡片EIGRL)
载荷集很容易地通过右键单击Model Browser窗口，并激活选项Create > Load Collector创建。 接着，赋予各自的名称与颜色。
图：创建SPC和力 (不需
Card image)
图：创建并命名 “buckling” 载荷集，使用Card Image EIGRL
点击Create将进入此载荷集的Card Editor
图：Card Image EIGRL
在此卡片中, 激活ND并输入1. 使RADIOSS 提取第一阶屈曲模态。
第三步
创建载荷和边界条件。确认各自的载荷集已设为当前（以SPC为例，在Model Browser中右键单击，并选择Make Current）。
图： “Model Browser”. 名字为SPC的载荷集设为当前并以粗体字显示
在约束面板中 (BCs > Create > Constraints) 选择梁底面的节点，如下图所示。约束所有自由度 (dof 1-6= 0)。
图：Constraint 面板。约束梁下端点并考虑其自由度情况。
接着，将力载荷集设为当前，交创建参考载荷（单位力）。在Forces面板中 (BCs > Create > Forces)选择梁的上端点，如下图示：
图： Force面板。施加一个沿X负向的参照单位力。(magnitude= -1, x-axis)
第四步
创建工况. 线性屈曲分析需要两个工况： 1. 线性静态加载 2. 线性屈曲分析
在Loadsteps面板中创建工况，可通过菜单栏并选择Setup > Create > LoadSteps进入。 确认type设为linear static。SPC在此调用单点约束载荷集（名字为SPC），且LOAD调用参照载荷（load）。
图：Loadsteps 面板(定义线性静态工况)
在第二个工况中, 屈曲工况,确认type设为linear buckling:
图：Loadstep面板(定义线性屈曲工况) STATSUB 调用静态分析结果且 METHOD调用屈曲模态信息的载何集
EIGRL
设置好了–到了运行分析的时候。
第五步
用RADIOSS开始计算。在Analysis页面，选择RADIOSS面板，选择工作目录并存储模型为.*fem (后缀名 *.fem指的是RADIOSS的输入文件)。最后，单击RADIOSS执行分析。
分析开始后，出现求解器视图。在Message log区域内显示计算步和可能出现的错误信息。
图：HyperWorks求解器页面
注：*.out文件包含计算前的模型检查信息和计算结果的其他基础信息。同时，如果出现错误可在此检查错误信息与警告，以助于调试输入文件。总之，即使分析成功完成也要经常地查看此文件。
第六步
分析后处理。如果计算成功，可看到工作目录下的新文件。 工作目录即模型存储目录。
打开HyperView查看结果。
确认并注意用以显示的工况。可在Results Browser中控制。
图：HyperView中的Subcase 2 “buckling”。蓝色为参考（未变形）模型
此杆计算出来的临界屈曲载荷(单位力1 N情况下) F = 863,53 N (临界屈曲因子λcr = 863.53)。换句话说，要发生屈曲，模型可承载863倍的单位载荷。当然，这只是个理想的线性屈曲分析，因其模型缺陷，真实的屈曲载荷一般会低于理想值。
例题: 机翼线性屈曲分析
此练习是对一个简化机翼进行线性屈曲分析。这在航空航天结构中是一个非常典型的问题，这类结构要求非常轻，因此就变得细薄。因为这些结构的细长比很高，屈曲失效验证成为必需。这个项目的目的是为了验证在对机翼施加3个静载荷工况下不会使其失效。
在此练习中，你将学会如何验证屈曲标准的基准设计。
模型信息 (Wing.hm)：
设计准则:
屈曲：1阶模态> (1.5 x).
(静态: U
铝材料:
ρ = 2.1e-9 T/mm3 [RHO]密度 E = 70.000 MPa [E] 杨氏模量 ν = 0.33 - [nu] 泊松比
在给定的模型中 (Wing.hm)，载荷，约束和工况已创建。在进行基准设计前，请先检查属性、材料和载荷集。 在此，以载荷集为重点：
图：左侧：在载荷集constraints中施加约束（绿色三角形）。右侧：在载荷集pressure中施加压强（黄色箭头）
图：翼尖载荷 (在载荷集tip_load中的紫色箭头)
在名为SUM的载荷集（其Card Image为LOAD）中，载荷集pressure (ID=2) 与tip_load (ID=3)中的载荷结合在一起。S，S1(1) ，和S1(2)是各自的比重因子。
三个静态工况为：
注：各分析的type都设为linear static。
为了定义屈曲工况，需要另一个载荷集。还记得那个遗忘的载荷集是哪个吗？
需要创建一个其Card Image为EIGRL的载荷集（定义感兴趣的屈曲模态的数目）。 其卡片信息为：
对于感兴趣的工况
下一步定义屈曲工况：
请比较所调用的载荷集（SPC= 和 METHOD =）和工况(STATSUB = )的ID
最后，运行分析并观察屈曲模态。
带重力载荷的屈曲
一般来讲，所有载荷相关的临界载荷因子在参照线性工况中定义（被STATSUB调用）。
假使临界载荷因子仅应用于可变载荷（例如压力、力等）而不能应用于重力载荷（恒定载荷），可以使用下面描述的两个方法、技术。下面描述的工作步骤见模型文件：Prestressed_gravity.hm.
工作方法A
假如考虑重力与可变载荷，需要以下载荷集：
• 可变载荷集（在此命名为force；无需Card Image） • 重力载荷集（Card Image为GRAV；需考虑系统单位）
•
将可变与恒定载荷结合在一起的载荷集（需要Card image为LOAD，在此命名为gravity_and_force）。
以上的卡片信息中，L1(1)链接重力载荷集，L1(2)指的是可变载荷（在这为force）。S,S1(1) 和S1(2)是各自的比例因子。
此外，要有一个模型的约束载荷集（命为spc，无Card image）和一个定义提取特征模态数目的载荷集（Card Image EIGRL，在这命名为EIGRL）。
第一步
确定所有载荷（恒定与可变载荷）的临界载荷因子。可将其称为buckling_all（即，运行考虑到所有载荷的屈曲分析）
第二步
调整恒定载荷因子，L1(1)，通过在载荷集的卡片信息中调整相关的比例值S1(1)。
S1(1) = 1 / [临界比例因子“buckling_all”]
可变载荷L1(2)的比例因子S1(2)保持不变，即S1(2)=1.0.
第三步
运行使用新的重力载荷L1(1)的比例值的屈曲分析。称这种额外屈曲分析和临界载荷因子为A_buckl_grav_force。
假如，屈曲因子A_buckl_grav_force等于buckling_all，一切好办了，因为如下：
重力载荷L1(1) = L1(1) / [临界载荷因子buckling_all] * [临界载荷因子A_buckl_grav_force]
然而，如果A_buckl_grav_force ≠ buckling_all 那么第二步与第三步必须不停重复。 现在gravity载荷L1(1)与A_buckl_grav_force成比例：
S1(1) = 1 / [A_buckl_grav_force]
这将导致另一个屈曲值，称之为B_buckl_grav_force，等等。
最终，当一个新的屈曲因子（例如，C_buckl …）与上一个（例如，B _buckl …）相等，一个额外用于检查的屈曲分析可完成。在此分析中，所有载荷与最后的屈曲因子，例如C_buckl …成比例。
此屈曲分析所描述的屈曲因子将接近于±1。
工作方法B
将RADIOSS求解器升级为HW11.0.220-HW (July 2012)，以上所描述的方法变得更简单直接。除载荷集spc与EIGRL (如上)外，现需要定义两个载荷集：
• 可变载荷集（在此名为force，无Card Image） • 重力载荷集(Card Image GRAV；注意系统单位)
基于这两个载荷集，定义两个工况： • 工况forces_only
• 工况gravity_only
注：现在，需要选择type: generic定义屈曲工况而不考虑重力、恒定载荷。 如前，SPC调用模型约束，METHOD (STRUCT)调用EIGRL。
STATSUB (BUCKLING) 调用forces_only，且
最后，重力载荷gravity_only以预载的形式施加。
可在prestressed_gravity.hm文件中查看载荷集与工况设定。