使用p型有限元法计算薄板结构屈曲载荷因子和相应屈曲形状的计算方法技术

本发明专利技术涉及一种使用p型有限元法计算薄板结构屈曲载荷因子和相应屈曲形状的计算方法，属于结构稳定性校核技术领域。本发明专利技术将p型有限元法应用于结构稳定性分析领域，求解在相应几何和约束条件下的应力场，得到相应的几何刚度矩阵并应用Lanczos迭代方法计算屈曲特征值，进而得到屈曲载荷因子和屈曲形状。该方法能减少计算成本，提高计算的收敛率和计算精度。

【技术实现步骤摘要】
使用p型有限元法计算薄板结构屈曲载荷因子和相应屈曲形状的计算方法
本专利技术涉及一种使用p型有限元法计算薄板结构屈曲载荷因子和相应屈曲形状的计算方法，属于结构稳定性校核

技术介绍
屈曲是结构常见的破坏形式之一，尤其是对于一些长、薄结构，具备这样性质的板壳结构广泛应用于航空、航天、航海、军工等领域中。飞机机翼、尾翼上的翼面壁板、梁腹板和机身上的蒙皮、隔框等常采用薄壁结构，这些结构都具有一个共同的特征，即一个方向的尺寸大大小于另外两个方向的尺寸，可以简化为板壳结构，当这种结构受到压缩、弯曲和剪切等外荷载单独或联合作用时，屈曲失稳是最常见的失效模式之一。薄板的屈曲分析在工程设计中至关重要，对于大多数薄板结构，屈曲失稳先于强度破坏。对于屈曲分析，主要应用伽辽金法、摄动法和有限元法。伽辽金法在求解过程中不涉及变分问题，简单易行，应用较广，但是在微分方程算子不是线性算子时很难求解；摄动法在求解过程中需要解渐近展开式，由于渐近级数一般是发散的，在求解时可能在级数展开项的前几项表现为逐渐收敛形式，而在级数项增加到一定程度后又表现为发散形式，要解决这一问题就需要较高的数学专业技能，增加了求解难度；有限元法采用将复杂结构离散化为一定数量的单元求数值解，为了提高数值解的计算精度，必须划分更小的单元。传统有限元法划分单元较多，计算效率偏低。本专利技术专利运用p型有限元法进行薄板结构的特征值屈曲分析，利用p型有限元法单元数量少、前处理少及收敛速度快的优点，提出了一种基于p型有限元法计算薄板结构屈曲载荷系数和相应屈曲形状的方法。
技术实现思路
本专利技术将p型有限元法应用于结构稳定性分析领域，求解在相应几何和约束条件下的应力场，得到相应的几何刚度矩阵并应用Lanczos迭代方法计算屈曲特征值，进而得到屈曲载荷因子和屈曲形状。该方法能减少计算成本，提高计算的收敛率和计算精度。本专利技术通过以下技术方案实现。使用p型有限元法计算薄板结构屈曲载荷因子和相应屈曲形状的计算方法，其步骤如下：步骤1、确定薄板结构的几何尺寸和材料参数，建立三维有限元模型，该步骤包括：步骤1.1、根据结构的实际受力情况，确定结构的受力特征；步骤1.2、根据结构所使用的材料确定材料属性，如有必要应通过实验测定材料的弹性模量、泊松比。步骤1.3、建立三维有限元模型。步骤1.4、施加载荷和边界条件：根据确定的载荷和边界条件，在相应的边界上施加相对应的载荷和边界条件；步骤2使用p型有限元法计算在对应荷载和约束条件下的线性解。具体步骤如下：步骤2.1、结合所指定的插值多项式阶数，依据方程Ka＝F(1)求解结构位移列阵a，其中，K为整体刚度矩阵，K＝∑GTKeG,Ke为单元刚度矩阵，a为结构位移列阵，F为结构节点载荷列阵，Ke，G，表示结构节点自由度与单元节点自由度的转换矩阵。步骤2.2、依据单元刚度矩阵公式单元等效节点载荷列阵其中：单元内部节点力外部节点力上式(2)至上式(5)中，Ω表示为单元内部；B＝LNI，B表示应变矩阵，L表示微分算子，NI表示高阶形函数矩阵，I表示单元的标号，Γt表示单元的外部边界；D表示应力矩阵，b表示体力，表示荷载边界条件，这些量由步骤1中确定的材料属性、载荷及边界条件计算得到；把上式(2)至上式(5)代入到上式(1)中，再对线性方程组(1)进行求解，则得到了结构位移列阵a；步骤2.3、依据公式aI＝Ga进一步得到了单元位移列阵aI；步骤2.4、依据公式u＝NIaI，ε＝Lu，σ＝Dε＝DBaI，求解得到应力场σ、应变场ε、位移场u。本专利技术采用p型有限元法求解位移场和应力应变场。p型有限元法可在固定网格下，通过采用阶谱类型的高阶形状函数来提高计算精度，其主要特征在于在三维标准六面体单元上采用了一种阶谱类型的高阶形状函数。本专利技术采用六面体网格单元对薄板屈曲载荷因子和屈曲形状进行数值模拟分析。阶谱类型的三维标准六面体单元高阶形状函数的构造形式如下：三维标准六面体单元如图1所示，其形函数为：(1)节点模式形状函数：节点模式形状函数共有8个，其定义如下：(2)边模式形状函数：边模式形状函数共有12(p-1)个，其定义如下：与连接节点1和节点2的边相关联的边模式形状函数为：其中，类似地，与连接节点2和节点3的边相关联的边模式形状函数为：其余边模式形状函数以此类推可得。(3)面模式形状函数：面模式形状函数共有3(p-2)(p-3),p≥4个，其定义如下：与连接节点1,2,5,6的面相关联的面模式形状函数为：其中，下标m＝m(i,j)由编号约定决定。其余面模式形状函数以此类推可得。(4)内部模式形状函数：内部模式形状函数共有p≥6个，其定义如下：其中，下标m＝m(i,j,k),i,j,k＝2,3,…,p-4,i+j+k＝6,7,…,p。插值多项式的阶次可依次从p＝1逐渐提高，提升插值多项式低阶的刚度矩阵可继续沿用，只需计算高阶部分，避免了刚度矩阵低阶部分的重复计算，具有良好的承袭性，节约了前处理的成本。步骤3：通过所求得的线性解中的应力场，计算所对应的几何刚度矩阵；再使用几何刚度矩阵进行屈曲载荷因子特征值求解：具体为步骤3.1、通过计算得到的应力场求得几何刚度矩阵；步骤3.2、将所求得的几何刚度矩阵用于特征值计算，将广义的特征值问题Kφ＝λMφ转化为标准特征值问题Kφ＝λφ，即将M化为单位阵；其中K为刚度矩阵，M质量矩阵，φ为特征矢量，λ即为所求得的屈曲载荷因子特征值；采用Lanczos迭代方法进行特征值计算；标准特征值问题：Kφ＝λφ任何初始向量U,||U||k＝1，U0＝0，其中||||k为k范数采用三项递推公式：{Uk+1}＝(K{Uk}-αk{Uk}-βk{Uk-1})/βk+1其中，β1＝0αk＝{Uk}TK{Uk}βk+1＝||K-{Uk}-αk{Uk}-βk{Uk-1}||2这里k＝1,2,3,...,m-1≤n，||||2代表2范数，m为最高阶特征值个数，n为方阵阶数，完成迭代计算后，最终得到两组解：用来建立预屈曲应力状态的线性解和特征值屈曲载荷因子解；步骤4：根据所求得的屈曲载荷因子是否满足精度要求，并判断能量范数误差是否满足一定的要求，如屈曲载荷因子未满足精度要求或者能量范数误差过大，提高插值多项式的阶数并增加网格剖分数量，返回步骤2；屈曲载荷因子满足精度要求，根据屈曲载荷因子获得相应屈曲形状。所述步骤2中p型有限元法在单元数量确定的情况下，通过提升插值多项式的阶数来提高计算结果的精度。当单元厚度较小时，使用低阶插值函数会使剪切刚度占主导地位，发生剪切锁定，因此插值多项式的阶本文档来自技高网...

【技术保护点】
1.一种使用p型有限元法计算薄板结构屈曲载荷因子和相应屈曲形状的计算方法，其特征在于步骤如下：/n步骤1：构建薄板构件三维有限元模型，该步首先创建薄板构件的几何模型，接着对所建立几何模型进行网格剖分；其次对要进行屈曲分析的薄板构件进行材料参数的设定；然后施加对应的荷载与边界条件；/n步骤2：使用p型有限元法计算在对应荷载和约束条件下的线性解中的应力场；/n步骤3：通过所求得的线性解中的应力场，计算所对应的几何刚度矩阵；再使用几何刚度矩阵进行屈曲载荷因子特征值求解：具体为/n步骤3.1、通过计算得到的应力场求得几何刚度矩阵；/n步骤3.2、将所求得的几何刚度矩阵用于特征值计算，将广义的特征值问题Kφ＝λMφ转化为标准特征值问题Kφ＝λφ，即将M化为单位阵；/n其中K为刚度矩阵，M质量矩阵，φ为特征矢量，λ即为所求得的屈曲载荷因子特征值；/n采用Lanczos迭代方法进行特征值计算；标准特征值问题：Kφ＝λφ/n任何初始向量U,||U||

【技术特征摘要】
1.一种使用p型有限元法计算薄板结构屈曲载荷因子和相应屈曲形状的计算方法，其特征在于步骤如下：
步骤1：构建薄板构件三维有限元模型，该步首先创建薄板构件的几何模型，接着对所建立几何模型进行网格剖分；其次对要进行屈曲分析的薄板构件进行材料参数的设定；然后施加对应的荷载与边界条件；
步骤2：使用p型有限元法计算在对应荷载和约束条件下的线性解中的应力场；
步骤3：通过所求得的线性解中的应力场，计算所对应的几何刚度矩阵；再使用几何刚度矩阵进行屈曲载荷因子特征值求解：具体为
步骤3.1、通过计算得到的应力场求得几何刚度矩阵；
步骤3.2、将所求得的几何刚度矩阵用于特征值计算，将广义的特征值问题Kφ＝λMφ转化为标准特征值问题Kφ＝λφ，即将M化为单位阵；
其中K为刚度矩阵，M质量矩阵，φ为特征矢量，λ即为所求得的屈曲载荷因子特征值；
采用Lanczos迭代方法进行特征值计算；标准特征值问题：Kφ＝λφ
任何初始向量U,||U||k＝1，U0＝0，其中||||k为k范数<...

【专利技术属性】
技术研发人员：张建铭，杨文升，陈峻，
申请(专利权)人：昆明理工大学，
类型：发明
国别省市：云南;53