一种基于提高迭代法稳定性和收敛性的位移求解方法.pdf

(19)中华人民共和国国家知识产权局(12)发明专利申请(10)申请公布号CN106126823A(43)申请公布日2016.11.16(21)申请号201610479906.5(22)申请日2016.06.23(71)申请人广州中国科学院工业技术研究院地址511458广东省广州市南沙区海滨路1121号(72)发明人丁桦刘红伟刘淑慧马荣鑫温靖张星(74)专利代理机构广州科粤专利商标代理有限公司44001代理人黄培智(51)Int.Cl.G06F17/50(2006.01)权利要求书1页说明书8页附图3页(54)发明名称一种基于提高迭代法稳定性和收敛性的位移求解方法(57)摘要本发明公开了一种基于提高迭代法稳定性和收敛性的位移求解方法，该方法基于有限元网格划分中包含奇异单元模型，从而单元刚度矩阵中会出现非对角占优的刚度矩阵为基础，对于非奇异单元及奇异单元中的非奇异节点采用传统的迭代技术进行处理，对于奇异单元中的畸形节点，可以在单元上进行改进，比如以预条件共轭梯度法为基础，采用体积作为惩罚函数进行预处理；本发明对非对角占优的有限元方程组进行修正，能够达到减少迭代次数，收敛速度更快，整个方法更加稳定的效果，其次采用先进行单刚和位移计算相乘后再叠加的方法，从而避免了组合总刚的计算，该计算量也是很庞大的；采用单刚的数据结构，方便数据的更新，并有利于并行计算。CN106126823ACN106126823A权利要求书1/1页1.一种基于提高迭代法稳定性和收敛性的位移求解方法，其特征在于，包括步骤：步骤1、取初始值u0,矩阵乘向量部分采用单元刚度矩阵乘以单元向量的做法，这样避免了合成总体刚度矩阵，则初始残余力向量r0＝F-∑eKeu0e上式中，F是载荷，Ke是单元刚度矩阵，u0e初始单元位移，r0是残余力向量；步骤2、解以下方程组Mz0＝r0求出z0，及后面迭代步中，同样解以下方程组求出zk；Mzk＝rk所述步骤2包括以下步骤：步骤21、对于非奇异单元及奇异单元中的非奇异节点，选取Jacobi预条件M＝diag(a11,a22,…ann)并记d＝M＝(d1,d2,…,dn)则，残差的求解如下：zk＝rk/dk这里k是第k次迭代，此处的预条件选取Jacobi-CG方法；步骤22、对于奇异节点，采用体积惩罚函数方法，体积惩罚函数求法如下：首先找到奇异节点，及其被哪些单元包含，求出该单元的质心，并找到奇异夹角所在的面，然后求出各个面的质心，质心坐标求法运用以下公式：从而求出这些面与对应单元质心组成该单元的体积，这就得到体积惩罚函数，记为v1e，v2e，…，vke，从而zki＝∑jvkjezkje这里zkje＝rkje/djei是奇异节点号，j是奇异节点所在单元的序列号；rkje是单元残差，dje是单元刚度矩阵对角元素。2.根据权利要求书1所述的基于提高迭代法稳定性和收敛性的位移求解方法，还包括步骤：步骤3、对k＝1，2，…计算y＝∑eKepkeuk＝uk-1+αkpkrk＝rk-1-αkyzki＝∑jvkjezkje步骤4、直到||rk||<ε，输出uk，结束。2CN106126823A说明书1/8页一种基于提高迭代法稳定性和收敛性的位移求解方法技术领域[0001]本发明涉及结构力学技术领域，具体涉及一种基于提高迭代法稳定性和收敛性的位移求解方法。背景技术[0002]通常在采用有限单元法对碰撞冲击问题进行数值分析时，特别是模拟高速冲击载荷下的结构响应时，例如弹打靶、爆炸过程中，一个实体变形前和拉伸变形后的情形如图6所示，节点A在进行弹打靶或爆炸前的位置如图1，一旦进行弹打靶或爆炸后，这个实体就发生了变形，继而点A的位置和方向就发生了变化，变成节点A’而线段就是节点A在变形过程中产生的位移，那么该节点的位移怎么求？[0003]如果要想描述整个物体位移，就要用有限单元的方法来描述，也即是用多个节点，及连接多个节点的网格(即是有限元网格)来代替原来的物体，那么求多个节点的位移就近似成求该物体的位移，也即是有限元的位移，那么该位移如何求？可以通过求解下面的有限元方程组来进行求解[0004]Ku＝F(1)[0005]其中K是单元刚度矩阵，F是物体在发生变形时所受到的力，u就是物体变形后的位置，即是位移。[0006]有限元网格往往会随着材料的剧烈变形而发生畸变，导致计算中出现非常畸形的网格如图3所示。或者对于比较复杂的结构在划分网格的过程中，就会形成非常畸形的网格，那么这个非常畸形的网格就是奇异的单元，从而方程组(1)的刚度矩阵就会出现非对角占优，甚至接近奇异的现象，从而方程组(1)就会出现条件数较差，求解困难的问题。那么对于包含这种非常畸形的网格的位移如何求解？现有的求解器方法如下：[0007]1)采用总体刚度矩阵的奇异值分解方法求出矩阵的奇异值，再进行求