在有限元分析中,为了计算各个单元的刚度矩阵和荷载向量,需要对单元上的积分进行数值近似。两种常见的积分方法是高斯积分和减缩积分。
高斯积分是一种数值积分方法,基于高斯-勒让德公式,用于近似计算函数在给定区间上的积分。在有限元分析中,高斯积分用于计算单元上的刚度矩阵和荷载向量。具体而言,对于每个单元,将其内部积分区域划分为一系列高斯积分点,然后根据高斯积分点的权重和位置,对积分表达式进行数值近似。通过对所有高斯积分点的贡献进行累加,可以得到单元上的刚度矩阵和荷载向量的近似值。
减缩积分(Reduced Integration)是一种在有限元分析中用于减少计算量和存储需求的技术,它采用了较少的高斯积分点来近似计算单元的刚度矩阵。减缩积分的核心思想是,通过在单元内部使用较少的高斯积分点来减少计算量,同时尽可能保持良好的数值稳定性。传统的高斯积分在每个单元上使用一组完整的高斯积分点,这可能会导致计算量过大,尤其是对于高阶元素或复杂几何形状的单元。减缩积分通过选择适当的积分点子集,例如在三维问题中选择积分点子集的体心点,来近似计算刚度矩阵。减缩积分的优点是减少了高斯积分点的数量,从而降低了计算量和存储需求。然而,它也可能导致计算结果的不精确性。由于减少了积分点的数量,减缩积分可能无法捕捉到一些细节信息,特别是在处理某些特殊情况(如奇异问题)时。总结来说,减缩积分是一种在有限元分析中用于减少计算量和存储需求的技术,它通过使用较少的高斯积分点来近似计算刚度矩阵。它可以在一定程度上降低计算开销,但也需要注意其对计算结果的影响。
1. 高斯积分与减缩积分的精度探究以平面应力问题为例,使用mathematical计算了四边形单元的刚度矩阵,对比了精确解,高斯积分解,减缩积分解的差异性。
Clear["Global`*"](*Parameters*)t=1;Em=1;v=0.2; (*mathematical计算默认都是精确解,但是我们提供的数据必须是分数;其次,对于大规模矩阵运算,精确解的求解计算非常慢*)Coord={{0,0},{3,2},{4,4},{2,5}};(*Constitutive Law-Plain Stress*)Cc=(Em/(1 v))*{{(1/(1-v)),(v/(1-v)),0},{(v/(1-v)),(1/(1-v)),0},{0,0,1/2}};(*Shape Functions*)N1=(1-eta)*(1-xi)/4;N2=(1-eta)*(1 xi)/4;N3=(1 eta)*(1 xi)/4;N4=(1 eta)*(1-xi)/4;Ni={N1,N2,N3,N4};Head[eta]; (*没有赋值直接拿来用的变量属于symbol数据类型*)Head[N4]; (*symbol数据类型的组合是Times数据类型,实际上还是symbol类型,Times代表symbol类型之间运算后赋值给N1*)(*Coordinate Mapping*)CoordMap={x,y};x=Sum[Ni[[i]]*Coord[[i,1]],{i,1,4}];y=Sum[Ni[[i]]*Coord[[i,2]],{i,1,4}];(*Jacobian Matrix*)NatCoord={xi,eta};J=Table[D[CoordMap[[i]],NatCoord[[j]]],{i,1,2},{j,1,2}];Jinv=Transpose[Inverse[J]];H=Jinv//MatrixForm ; (*这里H仅仅是矩阵格式,不是原来的Jinv,这种格式就是用来显示的,一种代表形式的数据类型*)(*B Matrix*)B1={{1,0,0,0},{0,0,0,1},{0,1,1,0}};B2={{Jinv[[1,1]],Jinv[[1,2]],0,0},{Jinv[[2,1]],Jinv[[2,2]],0,0},{0,0,Jinv[[1,1]],Jinv[[1,2]]},{0,0,Jinv[[2,1]],Jinv[[2,2]]}};B2//MatrixForm;B3=Table[0,{i,1,4},{j,1,8}]; (*如果是四边形八节点的二次单元,这里改为{j,1,8}*)Do[B3[[j,2i-1]]=B3[[j 2,2i]]=D[Ni[[i]],NatCoord[[j]]],{i,1,4},{j,1,2}];B3//MatrixForm; (*不要写成B=B1.B2.B3//MatrixForm*)B=B1.B2.B3;(*Stiffness Matrix*)IntFunc=(Transpose[B].Cc.B)*Det[J]; ke=t*NIntegrate[IntFunc,{xi,-1,1},{eta,-1,1}];ke//MatrixForm; (*如果计算精确解,用Integrate,即,kn=t*Integrate[IntFunc,{xi,-1,1},{eta,-1,1}],但是数据必须都是分数*)(*下面是用高斯积分,又算了一次刚度矩阵中的积分*)xiset={1/Sqrt[3],-1/Sqrt[3]};etaset={1/Sqrt[3],-1/Sqrt[3]};w={1,1}; (*权重系数列表*)kfull=t*Sum[w[[i]]*w[[j]]*(IntFunc/.{xi->xiset[[i]],eta->etaset[[j]]}),{i,1,2},{j,1,2}];(*下面是用减缩积分,又算了一次刚度矩阵中的积分*)kreduced=t*2*2 IntFunc/.{xi->0,eta->0};(*对比三种积分方法的差异*)ke//MatrixForm ;kfull//MatrixForm;kreduced//MatrixForm;
刚度矩阵精确结果
刚度矩阵高斯积分结果
刚度矩阵减缩积分结果
2. 结果分析与建议如上图所示,高斯积分的结果与解析解几乎一致,减缩积分的误差略大,可见减缩积分的确会降低精度。进一步改变单元形状发现,误差与单元不规则程度相关,例如,单元越接近矩形,精度越高。其中当单元形状为矩形时,高斯积分解等于精确解,这就是有限元软件中强调网格质量的原因(当然加密网格密度也是一种解决方法),或许很多人都因为单元扭曲严重而被终止或者警告分析。当然如果因为网格问题影响收敛时,不妨考虑以下方法:
重新网格划分:如果单元扭曲严重,可能需要重新划分网格,以改善单元的形状和质量。使用更合适的网格生成算法或软件工具,可以生成更好的网格布局,减少扭曲现象的发生。使用更适合的单元类型:对于特定类型的问题,某些单元类型可能更适合,并能够减少单元扭曲。局部修复或调整:对于少数扭曲严重的单元,可以尝试通过局部调整来修复它们。这可能包括重新定义节点位置或移动单元的某些节点,以改善其形状。增强积分(Enhanced Integration):增强积分是一种用于提高有限元分析精度的技术。通过在高应变或高梯度区域增加高斯积分点,可以提高对这些区域的数值近似。增强积分可以帮助解决扭曲导致的数值不稳定性问题。错误估计与自适应网格:使用错误估计技术来检测扭曲严重的单元,并基于估计的错误信息进行自适应网格细化或粗化。这可以帮助优化网格,减少扭曲并提高计算结果的准确性。3.结语最近刚复习了一遍有限元分析基础,旨在希望加深理解原理,进一步更好的应用有限元软件,有不当之处,欢迎指正,后续会抽时间更新其他笔记,感兴趣的同学可以看看我的上一篇有限元编程之杆单元,希望大家一起学习交流。