自编有限元程序中的刚度矩阵如何与Abaqus的刚度矩阵如何进行对照？

本次分享的是自编有限元程序中的单元刚度矩阵如何与Abaqus的单元刚度矩阵进行对照？

首先，讲一下为什么要进行单元刚度矩阵对比呢？

在自己的程序中，单元刚度矩阵的正确性直接影响到结果的精度，如果单元刚度矩阵没毛病，那组装也就是个局部自由度转换至全局自由度上的活儿，边界条件也是有固定的程序，就可以保证位移场的精度。这也就是UEL子程序的核心为单元刚度矩阵的计算的原因（自己瞎猜的！）。

再说为什么要和Abaqus进行对比呢？

因为木木是Abaqus的死忠粉儿！干啥都想对标Abaqus~

进入我们今天的正题！先和大家讲述一下Abaqus输出单元刚度矩阵的方法吧~

修改.inp文件

在inp语句*Output, field, variable=PRESELECT的下方加入：

*ELEMENT MATRIX OUTPUT, ELSET=element-50-91, STIFFNESS=YES, OUTPUT FILE=USER DEFINED, FILE NAME=ElementStiffness

即可输出单元集element-50-91的单元刚度矩阵，这个单元集 合可以自定义哈，这里我就是做了个演示~

然后会在工作目录生成ElementStiffness.mtx文件，可以用文本编辑器打开，

原以为是全矩阵的形式展现，实则：

矩阵的格式有点耐人寻味...

于是乎，今天的重点插件要来咯~

使用EMM插件

技术 邻的SnowWave02团队开发了一个小插件专门用于导出Abaqus中刚度矩阵，可以支持稀疏矩阵、全矩阵、单元刚度矩阵、整体刚度矩阵、质量矩阵、载荷矩阵，导出的形式为.mat文件，可以直接在Matlab中呈现，下载地址：https://www.jishulink.com/post/341364?toReply=true。里面包含用户帮助文档，很好操作的。

接下来试试自编程序与Abaqus内核求解的刚度矩阵对比，以三角形、四边形混合单元网格模型为例。左端固定，右端节点给一个正方形的位移，值为1，节点合位移对比图如下：

 

 

首先精度是可以保证的！然后我以2号单元和1009号单元为例，分别是四边形单元和三角形单元：

使用EMM插件将其刚度矩阵导出为.mat文件，在Matlab加载，可显示其刚度矩阵：

自编有限元程序的两个单元刚度矩阵显示为：

结果是一致的！

恭喜你，刚度矩阵组装之前的程序是正确的，如果程序的后处理结果还是不对，那就要检查一下边界条件咯~