【ANSYS APDL】提取结构刚度矩阵完整过程--附Python代码实现

Shihong Sun

已于 2022-10-25 00:22:10 修改

阅读量8.5k

点赞数 38

分类专栏： ANSYS APDL刚度矩阵文章标签： 1024程序员节矩阵

于 2022-10-24 22:30:46 首次发布

本文链接：https://blog.csdn.net/m0_61038364/article/details/127498110

版权

ANSYS APDL刚度矩阵专栏收录该内容

1 篇文章 1 订阅

订阅专栏

0 引言

结构刚度矩阵是什么

在这里简要回顾一下结构力学知识，可参考朱慈勉老师的《结构力学》

结构刚度矩阵来源于矩阵位移法，其中包括单元刚度矩阵，总刚度矩阵，结构刚度矩阵。在二维问题求解过程中，依次计算三者，最后求解可得到所有的节点位移U。该过程也称为直接刚度法，分析步骤如下：

第一步，首先介绍单元刚度方程的基本表达形式，单元刚度矩阵由材料力学表达，同理可加入转角，形成完整单元刚度矩阵。

第二步，将得到的所有单元刚度矩阵（注意需要经过局部系到整体系的转换）拼接，得到总刚度方程。

第三步，结合位移边界条件和外力，修正总刚度方程。位移约束用来修正U，若完全约束，则将对应行列删除。外力用来修正F，在相应节点的力定义为外力，其余力均为0。此时的刚度矩阵为非奇异，可以直接进行矩阵求解。

为什么要提取结构刚度矩阵

ANSYS用到的基本原理，就是上面提到的直接刚度法。在一些二维的问题中，将ANSYS的结构刚度矩阵提取出来，进行求解，得到的位移结果和ANSYS运行结果完全一样。因此采用刚度矩阵可以更直观地展示计算过程，提升运行效率。对于一些需要多次求解、变化约束、变化外力条件等情况，刚度矩阵的优势尤为明显。

1 ANSYS APDL如何导出结构刚度矩阵

导出步骤

这里直接给出导出刚度矩阵的命令（前提是前处理结束，划分完单元格，施加约束和外力）

!进行模态分析
/SOLU
/OUTPUT, 'Modal_res',dat
ANTYPE, MODAL
MODOPT, LANB, 20
MXPAND, 20, , ,NO
SOLVE
FINISH

接下来可采用HBMAT命令从“Model_1.full”中提取整体刚度矩阵

/AUX2
FILE,'Model_1',full
HBMAT, 'Stiffness_mat', dat, , ASCII, STIFF, YES, YES
HBMAT, 'Mass_mat', dat, ,ASCII, MASS, YES, YES
FINISH

其中，Model_1需要改为自己路径下的.full文件名。Stiffness代表生成刚度矩阵K，Mass代表生成结构整体质量矩阵。

文件解释

此时，我们的工作路径下生成了Stiffness_mat.dat文件和Stiffness_mat.mapping文件。文件可直接用记事本打开。

Stiffness_mat.mapping文件中列含义：第一列为方程编号，即带宽优化后矩阵的行/列编号；第二列和第三列分别指定与矩阵某一行/列对应的节点编号和自由度标签。

在这里插入图片描述
Stiffness_mat.dat文件含义。
第一行为文件内容的说明，指出该文件为ANSYS二进制文件导出的刚度矩阵，存储格式为Harwell-Boeing format。
第二行包含五个整数：第一个“98857”说明文件中具有98857个数据行；第二个数字“12961”说明数据的前12961行（文件的6-12966行）给定矩阵非零元素（对角线元素）的列索引指针（注意不是列索引），该数字等于自由度数加1；第三个整数“36468”指出文件中接下来的36468行给定非零元素的行索引；第四个整数“12960”指出文件接下来的12960行给出矩阵非零元素¹。

预处理

将Stiffness_mat.dat文件打开后，它的完整内容从上到下可以划分为三部分，这个很容易就区分出来了。我们要的是第二和第三部分。新建一个xlsx文件，命名为dat文件。按照下面图所示，把Stiffness_mat.dat的内容复制过来。
在这里插入图片描述

在运算之前检查一下最后面，数一下有几个数连续，记为count。

2制作结构刚度矩阵（python实现）

直接附python代码

import pandas as pd
import numpy as np
df=pd.read_excel(r'\dat文件.xlsx',sheet_name='dat文件',header=None)
def stiffness(p,q,count):
    '''p为`Stiffness_mat.dat`文件第二行第3个数，q为第二行最后一个数'''
    df_list=np.array(df.iloc[:p,0]).squeeze().tolist()
    result=[]
    n=1
    for i in range(1,p,1):
        if df_list[i]>df_list[i-1]:
            result.append(n)
        else:
            result.append(n)
            n=n+1
    result.append(0)
    for i in range(count):     
        result[-i-1]=df_list[-i-1]
            
    df_result=pd.DataFrame(np.zeros((q,q)))
    for i in range(p):
        x=df.iloc[i,0]
        y=result[i]
        df_result.iloc[int(x)-1,int(y)-1]=df.iloc[i,3]
    
    k_full=df_result+df_result.T
    for i in range(q):
        k_full.iloc[i,i]=df_result.iloc[i,i]
    return k_full

stiffness(36468,12960,6).to_csv(r'\结构刚度矩阵.csv',index=False,header=None)

Note

以上方法大多数情况通用。有个别特殊情况，需要对result手动编号。这个时候更改函数输出为result，手动调整最后几个特殊情况的数，再运行后续代码。
在这里插入图片描述

http://feacat.com/article1/ANSYS%20%E4%B8%AD%E6%95%B4%E4%BD%93%E5%88%9A%E5%BA%A6%E7%9F%A9%E9%98%B5%E7%9A%84%E8%BE%93%E5%87%BA%EF%BC%88%E4%B8%80%EF%BC%89%EF%BC%9A%E5%AD%90%E7%BB%93%E6%9E%84%E5%88%86%E6%9E%90.html ↩︎