最近在考慮自己編寫的程序和商用軟件的驗證問題，有限元結構分析中最關鍵的一環就是剛度矩陣的獲得，如果涉及到模態分析，還有質量矩陣。考慮到商業軟件的成熟性，可以用ANSYS生成的剛度矩陣做參照來看自己編寫的程序是否正確，因此如何提取ANSYS中結構的剛度矩陣，并進行隨后的驗證或者二次開發是一個問題。

http://www.yqgqt.org.cn/content/post/1796144

受上述帖子的啟發，使用MATLAB提取ANSYS中的剛度和質量矩陣，并進行模態分析驗證提取的矩陣的正確性。

首先，在ANSYS中使用HBMAT命令方法提取整體矩陣。

命令：HBMAT,fname,ext,--,form,matrx,rhs

其中：

Fname---輸出矩陣的路徑和文件名，缺省為當前工作路徑和當前工作文件名。

ext---輸出矩陣文件的擴展名，缺省為.matrix。

form---定義輸出矩陣文件的格式，其值可取：

=ASCII：ASCII碼格式；

=BIN：二進制格式。

matrix---定義輸出矩陣的類型，其值可取：

=STIFF：輸出剛度矩陣。可用于寫入了.FULL文件的任何類型的分析。

=MASS：輸出質量矩陣。可用于特征值屈曲、子結構分析、模態分析。

=DAMP：輸出阻尼矩陣。僅用于有阻尼的模態分析。

rhs---右邊項輸出控制（右邊項指用矩陣所表示方程的等號右端矢量，這里可為節點荷載向量），如rhs=YES則輸出，如rhs=NO則不輸出。

模態分析時，因僅LANB和QR法可生成完整的質量矩陣，因此也僅采用這兩種方法時才可使用HBMAT命令得到質量矩陣文件。

這個是Ansys幫助文件中對HBMAT命令的解釋，在Ansys分析中，會在根目錄下形成.full的二進制文件，里面存儲的就是已經組合好的剛度和質量矩陣，只是是以二進制文件的形式存在，通過HBMAT命令可以將二進制文件轉換為可以使用的十進制文件，十進制文件是以txt的形式存放在根目錄下。

ANSYS中根據Harwell-Boeing文件格式存儲剛度矩陣，如圖所示：

文件基本格式是前面有4或5行描述數據，其后為單列矩陣元素值，說明如下：

第1行：格式（A72），為文件頭的字符型解釋，如剛度矩陣或質量矩陣等標題。

第2行：格式（5I14），分別表示該文件的總行數（不包括文件頭）、矩陣列指針的總行數、矩陣行索引的總行數、矩陣元素數值的總行數、右邊項總行數。

第3行：格式（A3,11X,4I14），分別為矩陣類型、矩陣行數、矩陣列數、矩陣行索引數（對組裝后的矩陣，該值等于矩陣行索引數）、單元元素數（對組裝后的矩陣此值為0）。

第4行：格式（2A16,2A20），分別表示列指針格式、行索引格式、系數矩陣數值格式、右邊項數值格式。

第5行：格式（A3,11X,2I14），A3各列分別表示右邊項格式、應用高斯起始矢量、應用eXact求解矢量；兩個整數分別表示右邊項列數、行索引數。三個字符中的第1個字符可取：F---全部存貯（如節點荷載向量的全部元素）、M---與系數矩陣相同方法。

第6行后：矩陣元素值（單列）。

矩陣類型用3個字符表示，第1個字符可取：R---實數矩陣、C---復數矩陣、P---僅矩陣結構（無元素數值）；第2個字符可取：S---對稱矩陣、U---不對稱矩陣、H---Hermitian矩陣、Z---病態對稱矩陣；R---帶狀矩陣；第3個字符可取：A---組裝的矩陣、E---單元矩陣（未組裝）。對稱矩陣只存儲下三角元素，如結構剛度矩陣為對稱矩陣，Harwell-Boeing格式則僅記錄下三角元素。

通過上述解釋我們可以知道該剛度矩陣有1027個列指針，有7978個行索引，有1026個右邊項，矩陣的維度為1026×1026。

要將上述文件讀取到MATLAB中并轉化為常規的稀疏矩陣，需要了解Harwell-Boeing壓縮儲存格式及重構方法，這里面采用的是CSC方式，以列指針、行索引的方式去重構稀疏矩陣。

以三層框架結構為參考例子，其模態分析命令流如下：

fini

/clear

/filename,'Frame'    

/prep7

et,1,beam188

et,2,shell63

mp,ex,1,3.0e10

mp,nuxy,1,0.2

mp,dens,1,2500

r,1,0.15,0.15,0.15,0.15

sectype,1,beam,rect

secdata,0.8,0.8

k,1,

k,2,,,5

k,3,10,,5

k,4,10

k,5,,3.5

k,6,,3.5,5

k,7,10,3.5,5

k,8,10,3.5

k,9,5,3.5,,

k,10,5,3.5,5

k,11,,7

k,12,,7,5

k,13,10,7,5

k,14,10,7

k,15,5,7,,

k,16,5,7,5

k,17,,10.5

k,18,,10.5,5

k,19,10,10.5,5

k,20,10,10.5

k,21,5,10.5,,

k,22,5,10.5,5

l,1,5

l,2,6

l,3,7

l,4,8

l,5,11

l,6,12

l,7,13

l,8,14

l,11,17

l,12,18

l,13,19

l,14,20

a,5,6,10,9

a,7,8,9,10

a,11,12,16,15

a,13,14,15,16

a,17,18,22,21

a,19,20,21,22

lsel,all

latt,1,,1,,,,1

lesize,all,,,4

lmesh,all

asel,all

aatt,1,1,2

esize,,4

amesh,all

allsel

nummrg,all

numcmp,node

nsel,s,loc,y,0

d,all,all

fini

!模態分析

/solu

allsel

antype,modal

modopt,lanb,10

solve

finish

!(2). 模態分析完成后，在ansys中運行HBMAT命令得到十進制文件：

! 獲取矩陣十進制文件

/aux2

file,Frame,full !與文件名一致

hbmat,Stiff,txt,,ascii,stiffness,yes,yes   !提取剛度矩陣 存儲為Stiff.txt

hbmat,Mass,txt,,ascii,mass,yes,yes     !提取質量矩陣  存儲為Mass.txt

Finish

提取其前10階模態如下：

通過MATLAB對文件提取并處理后，可以使用spy命令查看稀疏矩陣的分布，剛度矩陣分布和各變量維度如下：

稀疏剛度矩陣的具體分布如下：

稀疏質量矩陣的具體分布如下：

計算矩陣的特征值，利用提取的剛度矩陣和質量矩陣做模態分析。

與ANSYS中的結果一致，因此驗證了正確性。

在Stiff.mapping文件中，格式是這樣的：

在MATLAB中使用相同的方法提取并整理，如下：