基于Python的SAP2000向ANSYS模型轉(zhuǎn)化技術及軟件開發(fā)研究

馬小平

(中鐵第一勘察設計院集團有限公司,西安 710043)

引言

隨著我國經(jīng)濟實力與科學技術的迅猛發(fā)展,一系列新穎復雜的建筑結構不斷涌現(xiàn),以滿足人民生活的需求。為保證此類復雜建筑結構的可靠性,工程人員常利用一些成熟的商業(yè)分析軟件對結構進行設計和校核。然而,結構設計是一個不斷反復調(diào)整以獲得最優(yōu)結果的過程,特別是對于復雜的大跨度或高層建筑,除常規(guī)的彈性分析之外,非線性分析也必不可少,例如:彈塑性時程分析、靜力穩(wěn)定分析、多尺度受力分析等[1-3]。因此,僅依靠單一的設計軟件或結構分析軟件,難以實現(xiàn)對復雜結構多類型力學響應的綜合分析。

工程人員為驗證設計結構的正確性和可靠性,常采用多種類型的軟件進行比較分析。除常規(guī)的設計軟件(如:PKPM和YJK)之外,SAP2000、ETABS等結構專業(yè)設計分析軟件,因其良好的三維結構分析能力、便捷的空間建模能力以及完善的荷載計算功能,受到了業(yè)內(nèi)專業(yè)人士的廣泛應用[4-6]。同時,ANSYS、ABAQUS等大型通用有限元分析軟件,因其超強的拓展性和靈活性,加上內(nèi)嵌精確的數(shù)值計算方法,故在特殊構件、復雜結構的分析驗算中得到使用[7-10]。然而,僅依靠ANSYS、ABAQUS自帶的前處理功能對較為復雜的結構模型進行重構,不僅耗時費力,且意義不大。考慮各類型軟件模型建立方法上的相似性,故可編譯軟件之間數(shù)據(jù)轉(zhuǎn)換的接口程序,將人工建模的工作量轉(zhuǎn)移到程序完成,極大地提高建模效率和準確性[11-13]。

結構設計分析軟件與通用有限元軟件之間并不能直接進行數(shù)據(jù)交換共享,因此,眾多企業(yè)和科研團隊基于各類開發(fā)平臺推出了一系列模型轉(zhuǎn)化程序。北京盈建科軟件股份有限公司已相繼開發(fā)了YJK和ETABS接口軟件(YJK-ETABS)、YJK和MIDAS接口軟件(YJK-MIDAS)、YJK和SAP2000接口軟件(YJK-SAP2000)、YJK和ABAQUS接口軟件(YJK-ABAQUS)等,為YJK設計模型向各類型結構分析軟件的轉(zhuǎn)化提供了便捷。葛金剛等[14]基于Python語言,開發(fā)了將結構模型由SAP2000及MIDAS向ABAQUS轉(zhuǎn)換的接口程序,并以天津某展覽館屋蓋結構為例,對比驗證程序轉(zhuǎn)化結果的正確性。孟仲永[15]以AutoCAD作為圖形處理平臺,利用ObjectArx及VS2008和C++語言開發(fā)了從SAP2000到ABAQUS的模型轉(zhuǎn)換程序,并應用到了兩個實際大跨度復雜空間結構的模型轉(zhuǎn)化中,驗證了轉(zhuǎn)化結果可靠性。在此基礎上,王杰[16]通過類似的方法研發(fā)了MIDAS/GEN到ABAQUS的模型轉(zhuǎn)換程序,實現(xiàn)了結構樓層信息、節(jié)點及單元信息、截面及材料信息、荷載及約束信息的高效轉(zhuǎn)化?；赟AP2000和ABAQUS軟件之間的內(nèi)在轉(zhuǎn)換邏輯,張月強等[17]編制了用以模型自動轉(zhuǎn)換的MTR1.0程序,實現(xiàn)了SAP2000幾何模型向ABAQUS中的精確轉(zhuǎn)化。基于此,祝輝慶等[18]也利用MATLAB軟件編制SAP2000到ABAQUS的模型數(shù)據(jù)轉(zhuǎn)換程序,實現(xiàn)了節(jié)點、單元、材料、截面、荷載、約束等信息的自動轉(zhuǎn)化。曹偉良等[19]采用VB6.0結合ACCESS數(shù)據(jù)庫程序技術,將SAP2000導出的MDB文件轉(zhuǎn)化為ANSYS命令流文件,實現(xiàn)了SAP2000常規(guī)構件及簡單荷載類型的轉(zhuǎn)化。

通過對現(xiàn)有模型轉(zhuǎn)化技術及軟件開發(fā)研究的深入調(diào)研,可以看出,眾學者分享了不同軟件之間的轉(zhuǎn)化思想以及相應程序平臺的開發(fā),為相關研究的進一步發(fā)展提供了有力的技術支撐。然而,針對結構設計分析軟件SAP2000與通用有限元軟件ANSYS之間的模型轉(zhuǎn)化,其相關研究整體上仍較為匱乏,且存在以下重要問題亟待改進。

(1)現(xiàn)有SAP2000向ANSYS的模型轉(zhuǎn)化程序,其轉(zhuǎn)化內(nèi)容過于簡單。部分自編程序僅支持幾何模型的轉(zhuǎn)化,相對成熟的程序雖然轉(zhuǎn)化內(nèi)容較為全面,但也僅支持簡單的低階單元(如Beam44梁單元、Shell63殼單元等)、常規(guī)截面類型(如矩形、工字形等)、常規(guī)荷載(均布荷載、集中荷載等)。當需要依賴高階單元(Beam188、Shell181等)確保非線性分析的準確性時,或結構包含特殊截面(如變截面、自定義截面等)和荷載類型(三角形分布荷載、梯形荷載等),上述轉(zhuǎn)化程序就會失效。

(2)現(xiàn)有SAP2000向ANSYS的模型轉(zhuǎn)化程序,其功能性嚴重不足。對于實際結構分析中常用到連接件設置、節(jié)點坐標系變換、梁單元坐標系變換、梁端自由度釋放、剛性域等功能,現(xiàn)有轉(zhuǎn)化程序仍然不夠完善。此外,由于未考慮控制網(wǎng)格劃分數(shù)量的功能,對計算精度或效率有要求的模型,現(xiàn)有轉(zhuǎn)化程序難以提供有效幫助。

(3)現(xiàn)有SAP2000向ANSYS的模型轉(zhuǎn)化程序,其設計思想過于簡單?，F(xiàn)有轉(zhuǎn)化程序大多以ANSYS中的節(jié)點(Node)代替SAP2000中結點(Joint),直接建立構件或結構的有限元模型。該方法雖然能夠快速實現(xiàn)模型的創(chuàng)建,但對于有網(wǎng)格劃分要求或荷載施加形式復雜的模型,以上述思想設計的程序無法實現(xiàn)模型的準確轉(zhuǎn)化。

鑒于此,開發(fā)團隊確定了以關鍵點作為模型轉(zhuǎn)化基礎,并考慮了多種單元類型、材料本構、荷載形式的轉(zhuǎn)化需求,研究了SAP2000模型結構信息在ANSYS中的實現(xiàn)條件及方法,進而開發(fā)了“SAP2000 To ANSYS模型轉(zhuǎn)化軟件 V1.0”(簡稱“STAMT V1.0”),以期實現(xiàn)SAP2000結構分析模型向ANSYS有限元仿真模型的全方面轉(zhuǎn)化,為結構設計的從業(yè)人員提供必要的技術支撐。

1 SAP2000 To ANSYS模型轉(zhuǎn)化程序

1.1 STAMT程序的目標及功能

STAMT程序的核心目標是依據(jù)SAP2000導出的.s2k文件,自動讀取文件內(nèi)數(shù)據(jù),按照指定的網(wǎng)格劃分數(shù),轉(zhuǎn)化模型的全部信息,生成ANSYS有限元軟件可讀取的參數(shù)化命令流,并以文檔的形式導出,實現(xiàn)SAP2000結構分析模型向ANSYS有限元分析模型的精準轉(zhuǎn)化。

STAMT程序具備了以下幾個主要功能。

(1)SAP2000文件的導入和數(shù)據(jù)讀取功能,如:材料信息讀取、幾何坐標讀取、構件截面讀取、荷載分布讀取等。

(2)自定義框架梁柱單元的網(wǎng)格劃分數(shù)量,即用戶可根據(jù)需求指定梁單元網(wǎng)格劃分的數(shù)量。

(3)多種結構單元類型及基本材料屬性的轉(zhuǎn)化,如:梁單元、板殼單元、彈簧單元、質(zhì)量單元、線彈性材料、彈塑性材料等。

(4)多種結構構件功能性需求的轉(zhuǎn)化,如:節(jié)點坐標系變換、梁單元坐標系變換、梁端自由度釋放、創(chuàng)建剛性域等。

(5)多種荷載工況的轉(zhuǎn)化,如:節(jié)點力、節(jié)點位移、梁柱分布載荷、梁柱集中載荷、板殼均布載荷、傳遞到梁的板殼載荷等。

STAMT程序面向的用戶主要是土木工程領域相關的結構設計及分析人員。因此,程序應具有良好的用戶交互界面,方便不熟悉ANSYS參數(shù)化設計語言(ANSYS Parametric Design Language,APDL)的人員進行系統(tǒng)應用,且轉(zhuǎn)化結果清晰直觀,可供相關人員快速、高效地實現(xiàn)模型的轉(zhuǎn)化,節(jié)省模型重建所耗費的時間,大大提高工程人員的結構驗核效率。

1.2 STAMT程序的開發(fā)環(huán)境

為滿足上述目標及功能,STAMT程序需具備良好的用戶交互界面,快速實現(xiàn)文件的導入導出,并具備模型數(shù)據(jù)的批量轉(zhuǎn)化和數(shù)值計算功能。Python作為一款解釋性、編譯性、面向?qū)ο蟮哪_本語言,不僅具有易于學習、閱讀、維護等特點,且包含了各類的標準庫,可實現(xiàn)不同文本數(shù)據(jù)的讀入和寫出,以及各類函數(shù)的數(shù)值運算[20-21]。同時,Python中的PyQt5庫涵蓋了豐富的功能函數(shù)用于交互界面設計,借助Eric6和Qt Designer可更加快速、便捷地實現(xiàn)軟件操作界面的創(chuàng)建。因此,本軟件采用Python 3.9作為開發(fā)環(huán)境,借助Eric6集成開發(fā)軟件,采用用戶界面與業(yè)務邏輯分離思想進行軟件的整體架構設計。

1.3 STAMT程序的關鍵技術及流程

STAMT程序內(nèi)部的模型轉(zhuǎn)化過程主要包括讀取SAP2000文件,生成并導出相關的APDL代碼,如定義單元類型、材料屬性、實常數(shù)、截面屬性等,創(chuàng)建關鍵點、幾何對象、單元等,變換單元坐標軸及施加約束,施加荷載及質(zhì)量等。具體轉(zhuǎn)化流程如圖1所示。

圖1 STAMT程序內(nèi)部模型轉(zhuǎn)化流程Fig.1 Model transformation process of the STAMT program

1.3.1 讀取SAP2000文件

SAP2000軟件可導出.s2k格式文件,該文件涵蓋了SAP2000模型的所有數(shù)據(jù)庫表格,包含了所有在交互界面中設置的模型信息。STAMT程序會首先讀取該文件,將數(shù)據(jù)分類并以列表或字典格式存儲。

1.3.2 單元類型的轉(zhuǎn)化

為滿足更高的計算需求,STAMT程序提供了梁單元(Beam188)、殼單元(Shell181)、彈簧單元(Combin14及Combin39)、質(zhì)量單元(Mass21)、網(wǎng)格單元(Mesh200)的轉(zhuǎn)化條件。其中,殼單元通過關鍵項KEYOPT(1)分為兩類,即考慮薄膜及彎曲剛度的殼單元,其對應SAP2000中的殼單元,以及僅考慮薄膜剛度的殼單元,其對應SAP2000中的膜單元。彈簧單元包含線性彈簧單元(Combin14)和非線性彈簧單元(Combin39)。通過Combin14的關鍵項KEYOPT(2),分別定義了6個自由度方向上的線性彈簧單元。通過Combin39的關鍵項KEYOPT(1)和KEYOPT(2),在6個自由度方向上分別定義了彈性非線性彈簧單元以及彈塑性非線性彈簧單元。以彈簧單元為例,首先,建立轉(zhuǎn)化彈簧單元的子函數(shù)。然后,依次編寫一維線性彈簧、一維彈性非線性彈簧、一維彈塑性非線性彈簧的轉(zhuǎn)化代碼。最后,將列表返回并調(diào)用子函數(shù)。

def Spr_Elem_Type():

APDL_Spr_Elem_Type+=['et,2,combin14'+' '+'keyopt,2,1,0'+' '+'keyopt,2,2,1'+' ']

APDL_Spr_Elem_Type+=['et,8,combin39'+' '+'keyopt,8,1,0'+' '+'keyopt,8,2,0'+' '+'keyopt,8,3,1'+' ']

APDL_Spr_Elem_Type+=['et,14,combin39'+' '+'keyopt,14,1,1'+' '+'keyopt,14,2,0'+' '+'keyopt,14,3,1'+' ']

return APDL_Spr_Elem_Type

1.3.3 材料屬性的轉(zhuǎn)化

為使轉(zhuǎn)化后模型的本構關系更具通用性,STAMT程序提供了多種基本材料屬性的轉(zhuǎn)化條件,包括:彈性、理想彈塑性、雙線性彈塑性、多線性彈塑性等。

1.3.4 實常數(shù)的轉(zhuǎn)化

STAMT程序提供了梁單元、殼單元、彈簧單元、質(zhì)量單元的實常數(shù)轉(zhuǎn)化條件。梁單元、殼單元僅創(chuàng)建了實常數(shù)編號,以滿足網(wǎng)格劃分的需要,無其他實際意義。彈簧單元的實常數(shù)包含了線性彈簧單元實常數(shù)、非線性彈性彈簧單元實常數(shù)、非線性塑性彈簧單元實常數(shù)三部分,具體設置參數(shù)和對應關系如表1所示。

表1 彈簧單元的對應關系及實常數(shù)設置Table 1 Correspondence of spring units and real constants

質(zhì)量單元的實常數(shù)可設置3個平動方向的質(zhì)量以及3個轉(zhuǎn)動方向的質(zhì)量慣性矩。確定了模型信息中實常數(shù)的對應關系,以及整體的轉(zhuǎn)化邏輯,可編譯相關代碼實現(xiàn)實常數(shù)信息的轉(zhuǎn)化。

1.3.5 截面屬性的轉(zhuǎn)化

STAMT程序提供了梁殼單元截面屬性定義的轉(zhuǎn)化條件。殼單元截面僅需設置板殼的厚度。梁單元截面可轉(zhuǎn)化的類型有矩形、工字形、箱形、圓管、自定義截面和變截面。其中,自定義截面和變截面的轉(zhuǎn)換較為特殊,現(xiàn)有轉(zhuǎn)化程序都未涉及,而這兩種類型在實際工程中應用廣泛。為此,結合ANSYS中有關自定義梁截面和變截面梁的APDL命令,編譯可轉(zhuǎn)化此類截面的相關代碼。

ANSYS中自定義截面的APDL命令如下:

APDL_Beam_Sec+=['sectype,'+Frame_Sec_Num+',beam,asec,'+Frame_Sec_Name+' '+'secdata,'+ASEC_Area[1]+','+ASEC_Iyy[1]+','+ASEC_Iyz[1]+','+ASEC_Izz[1]+','+ASEC_Iw+','+ASEC_Jt[1]+','+ASEC_CGy+','+ASEC_CGz+','+ASEC_SHy+','+ASEC_SHz+','+ASEC_TKz+','+ASEC_TKy+','+ASEC_TSxz+','+ASEC_TSxy]

上述代碼中,Frame_Sec_Num為截面編號;Frame_Sec_Name為截面名稱;ASEC_Area為自定義截面面積;ASEC_Iyy為自定義截面對y軸的慣性矩;ASEC_Iyz為自定義截面的慣性積;ASEC_Izz為自定義截面對z軸的慣性矩;ASEC_Jt為自定義截面的扭轉(zhuǎn)慣性矩。上述參數(shù)在SAP2000導出的數(shù)據(jù)文件中均有對應信息,可直接轉(zhuǎn)化得到。然而,自定義截面的翹曲慣性矩ASEC_Iw、截面重心y坐標ASEC_CGy、截面重心z坐標ASEC_CGz、截面剪切中心y坐標ASEC_SHy、截面剪切中心z坐標ASEC_SHz、截面沿z軸厚度ASEC_TKz、截面沿y軸厚度ASEC_TKy、截面xz剪切修正系數(shù)ASEC_TSxz、截面xy剪切修正系數(shù)ASEC_TSxy,無法從SAP2000文件中直接獲取并轉(zhuǎn)化得到。因此,程序在輸出的文件中保留有參數(shù)符號及說明,使用者可利用其他截面計算軟件求得上述參數(shù)后,填寫入命令流中進行仿真計算。

ANSYS中變截面的APDL命令如下:

APDL_Frame_Sec+=['sectype,'+Frame_Sec_Num+',taper,,'+Frame_Sec_Name+' '+'secdata,'+Sec_Num_Start+','+Joint_I_X+','+Joint_I_Y+','+Joint_I_Z+' '+'secdata,'+Sec_Num_End+','+Joint_J_X+','+Joint_J_Y+','+Joint_J_Z]

如上述代碼所示,首先,需要確定變截面梁的I節(jié)點(Frame_Joint_I)和J節(jié)點(Frame_Joint_J)的編號,并根據(jù)節(jié)點坐標數(shù)據(jù)(Joint_Coor)分別遍歷得到I、J節(jié)點的坐標。然后,在截面數(shù)據(jù)中,確定梁始端的截面編號(Sec_Num_Start)及其對應I節(jié)點坐標,梁末端的截面編號(Sec_Num_End)及其對應J節(jié)點坐標。通過以上步驟,即可轉(zhuǎn)化得到定義梁變截面屬性的APDL命令流。

1.3.6 關鍵點、幾何對象、單元劃分的轉(zhuǎn)化

利用SAP2000文件中的“連接數(shù)據(jù)”,通過編譯轉(zhuǎn)化代碼,可實現(xiàn)幾何模型(包括點、線、面)的創(chuàng)建。此處代碼無特殊難點,故不再贅述。STAMT程序提供了梁單元劃分、殼單元劃分、彈簧單元創(chuàng)建的轉(zhuǎn)化能力。其中,梁單元劃分的數(shù)量可根據(jù)需要輸入,殼單元按默認方式劃分。彈簧單元創(chuàng)建的代碼較為復雜,此處將詳細闡述。

首先,在連接(Link)位置上創(chuàng)建幾何線,選擇線并賦予空的材料屬性和實常數(shù),以及網(wǎng)格單元類型Mesh200,并劃分一個單元。

APDL_Spr_Elem+=['l,'+Link_Joint_I+','+Link_Joint_J]

APDL_Spr_Elem+=['lsel,s,line,,Line_Num_Max+1']

APDL_Spr_Elem+=['latt,Mat_Num_Max+1,Real_Num_Max+1,23']

APDL_Spr_Elem+=['lesize,all,,,1']

APDL_Spr_Elem+=['lmesh,all']

以x方向一維線性彈簧為例,確定單元類型和實常數(shù)后,即可根據(jù)連接(Link)節(jié)點I和J的坐標創(chuàng)建彈簧單元。其他單元類似。

APDL_Spr_Elem+=['type,2']

APDL_Spr_Elem+=['real,'+Real_Link_Linear]

APDL_Spr_Elem+=['en,Elem_Num_Max+1,node(kx('+Link_Joint_I+')'+','+'ky('+Link_Joint_I+')'+','+'kz('+Link_Joint_I+')'+')'+',node('+'kx('+Link_Joint_J+')'+','+'ky('+Link_Joint_J+')'+','+'kz('+Link_Joint_J+')'+')'+' ']

1.3.7 坐標軸變化及約束的轉(zhuǎn)化

STAMT程序提供了節(jié)點坐標軸變化、梁端自由度釋放、梁單元坐標軸變化、創(chuàng)建剛性域、施加支座約束等功能的轉(zhuǎn)化條件。依據(jù)SAP2000文件中“節(jié)點局部軸指定1-標準”“框架釋放指定1-通用”“節(jié)點約束指定”“節(jié)點剛性約束指定”“框架局部軸指定1-標準”等模塊的數(shù)據(jù),即可編譯代碼實現(xiàn)相應模型信息的轉(zhuǎn)化。節(jié)點坐標軸變換、梁端自由度釋放、施加支座約束分別通過APDL中的nmodif命令、endrelease命令、d命令等實現(xiàn),其代碼較為簡單,此處不再贅述。此處,重點介紹梁單元坐標軸變化及剛性域創(chuàng)建。

梁單元坐標軸變化的部分代碼如下。

首先,確定梁單元坐標系x軸的方向向量,并計算向量的模和歸一化后的單位向量。

Frame_x_Axis_Vec=np.array([(Joint_J_X-Joint_I_X),(Joint_J_Y-Joint_I_Y),(Joint_J_Z-Joint_I_Z)])

Frame_x_Axis_Vec_Len=np.linalg.norm(Frame_x_Axis_Vec)

Frame_x_Axis_Vec_Norm=Frame_x_Axis_Vec/Frame_x_Axis_Vec_Len

然后,將x軸的單位方向向量與整體坐標系下Z軸的方向向量作比較,若二者平行,則單元坐標系的y軸方向與整體坐標系下Y軸的方向一致,若不平行,則通過二者叉積確定。單元坐標系z軸的方向則通過單元坐標系的y軸單位方向向量與x軸單位方向向量的叉積得到。

Global_Z_Axis=np.array([0,0,1])

if Frame_x_Axis_Vec_Norm.dot(Global_Z_Axis)/(np.sqrt(Frame_x_Axis_Vec_Norm.dot(Frame_x_Axis_Vec_Norm))*np.sqrt(Global_Z_Axis.dot(Global_Z_Axis)))==1 or Frame_x_Axis_Vec_Norm.dot(Global_Z_Axis)/(np.sqrt(Frame_x_Axis_Vec_Norm.dot(Frame_x_Axis_Vec_Norm))*np.sqrt(Global_Z_Axis.dot(Global_Z_Axis)))==-1:

Frame_y_Axis_Vec=np.array([0,1,0])

Frame_z_Axis_Vec=np.cross(Frame_x_Axis_Vec_Norm,Frame_y_Axis_Vec)

else:

Frame_y_Axis_Vec=np.cross(Global_Z_Axis,Frame_x_Axis_Vec_Norm)

Frame_z_Axis_Vec=np.cross(Frame_x_Axis_Vec_Norm,Frame_y_Axis_Vec)

接著,確定旋轉(zhuǎn)中心軸(Rot_Axis)為單元坐標軸x,旋轉(zhuǎn)角度為Rot_Ang,起始向量(Vec_Origin)為單元坐標軸z,則可計算得到目標向量(Vec_Target),進而將目標向量的元素與梁單元起始節(jié)點(I節(jié)點)的坐標對應相加,即可得到旋轉(zhuǎn)后梁單元的方向節(jié)點。

Rot_Axis=Frame_x_Axis_Vec_Norm

Vec_Origin=Frame_z_Axis_Vec

Vec_Target=Vec_Origin*math.cos(Rot_Ang)+np.cross(Rot_Axis,Vec_Origin)*math.sin(Rot_Ang)+Rot_Axis*(np.dot(Rot_Axis,Vec_Origin))*(1-math.cos(Rot_Ang))

New_Frame_Z_Axis_X_Coor=Vec_Target[0]+Joint_I_X

New_Frame_Z_Axis_Y_Coor=Vec_Target[1]+Joint_I_Y

New_Frame_Z_Axis_Z_Coor=Vec_Target[2]+Joint_I_Z

最后,通過創(chuàng)建新的方向控制節(jié)點,并利用emodif命令修改梁單元的屬性,實現(xiàn)梁單元坐標軸的變化。

APDL_Beam_Axis+=['allsel,all']

APDL_Beam_Axis+=['n,'+'Node_Num_Max+'+str(i+1)+','+str(New_Frame_Z_Axis_X_Coor)+','+str(New_Frame_Z_Axis_Y_Coor)+','+str(New_Frame_Z_Axis_Z_Coor)]

APDL_Beam_Axis+=['lsel,s,line,,'+Frame_Line_Name]

APDL_Beam_Axis+=['esll,s']

APDL_Beam_Axis+=['emodif,all,-3,Node_Num_Max+'+str(i+1)+' ']

剛性域創(chuàng)建的部分代碼如下。首先,根據(jù)節(jié)點約束信息(Constraint_Joint),明確需要創(chuàng)建剛性體約束的節(jié)點,將其篩選出來并定義為一個組件(Node_Rigid_Region)。

for k in range(0,len(Constraint_Joint)):

if k==0:

APDL_Rigid_Region+=['ksel,s,kp,,'+Constraint_Joint[k]]

elif k>0:

APDL_Rigid_Region+=['ksel,a,kp,,'+Constraint_Joint[k]]

APDL_Rigid_Region+=['nslk,s']

APDL_Rigid_Region+=['cm,Node_Rigid_Region'+str(i+1)+',node']

接著,從每個剛性體約束節(jié)點組中選第一個節(jié)點作為主節(jié)點(Node_Main),并利用cerig命令在主結構與其他節(jié)點之間建立剛性約束。

APDL_Rigid_Region+=['Node_Main='+'node('+'kx('+Constraint_Joint[0]+'),ky('+Constraint_Joint[0]+'),kz('+Constraint_Joint[0]+')'+')']

APDL_Rigid_Region+=['cmsel,s,Node_Rigid_Region'+str(i+1)+',node']

APDL_Rigid_Region+=['cerig,Node_Main,all,'+'ux']

1.3.8 荷載及質(zhì)量的轉(zhuǎn)化

STAMT程序提供了節(jié)點質(zhì)量、節(jié)點荷載、支座位移、梁柱分布荷載、梁柱集中荷載、樓面均布荷載、樓面?zhèn)鬟f至梁的荷載等轉(zhuǎn)化條件。節(jié)點質(zhì)量、節(jié)點荷載、支座位移的轉(zhuǎn)化較為簡單,此處不再詳述。下面以梁柱荷載和樓面荷載的轉(zhuǎn)化作為關鍵技術點,進行展開介紹。

以梁柱分布荷載的轉(zhuǎn)化為例,由于對模型進行了網(wǎng)格劃分,因此框架梁上的分布荷載并不能直接通過sfbeam命令進行施加。鑒于此,分4個步驟編譯了梁單元分布荷載施加的轉(zhuǎn)化代碼。第一步,選取需要施加分布荷載的框架梁(Frame),確定其起始節(jié)點的幾何坐標以及梁長。第二步,基于Python語言,編譯APDL命令流,命令流內(nèi)容為選取ANSYS模型中梁(Beam)上所有節(jié)點,將其錄入至所創(chuàng)建列表中。第三步,與第二步方法一致,在ANSYS中遍歷每個梁上的節(jié)點,并將其定義為節(jié)點組件。以上步驟轉(zhuǎn)化代碼不再詳述。第四步,讀取分布荷載數(shù)據(jù),判斷荷載方向,確定其與單元坐標軸的夾角,最后施加于梁單元上。此處,著重闡述第四步的實現(xiàn)方法。

以單元坐標軸x和荷載沿重力方向為例,其他方向類似。首先,判斷框架分布荷載的方向,定義荷載的方向向量。

if Frame_Distr_Load_Dir=='Gravity’ or Frame_Distr_Load_Dir=='Gravity Proiected':

Load_Vec=np.array([0,0,-1])

然后,判斷梁單元坐標方向與荷載方向的夾角。單元坐標軸的方向向量(Frame_x_Axis_Vec)已在前文說明。

Frame_x_Axis_Vec_Len=np.sqrt(Frame_x_Axis_Vec.dot(Frame_x_Axis_Vec))

Cos_Angle_x=Frame_x_Axis_Vec.dot(Load_Vec)/(Frame_x_Axis_Vec_Len*Load_Vec_Len)

Angle_x=np.arccos(Cos_Angle_x)*180/np.pi

進而,判斷框架分布荷載是否為投影荷載,若為投影荷載,則根據(jù)SAP2000程序的規(guī)定,將框架分布荷載乘以荷載與坐標軸夾角的正弦值。

if Frame_Distr_Load_Dir=='Gravity Proiected':

Frame_Distr_Load_A=str(float(Frame_Distr_Load_A)*math.sin(Angle_x))

Frame_Distr_Load_B=str(float(Frame_Distr_Load_B)*math.sin(Angle_x))

最后,將框架分布荷載在單元坐標軸x上投影,并根據(jù)分布荷載在框架梁上的作用距離,篩選梁上節(jié)點,分段將荷載施加到梁單元上。

Load_a_x=float(Frame_Distr_Load_A)*math.cos(Angle_x)

Load_b_x=float(Frame_Distr_Load_B)*math.cos(Angle_x)

APDL_Beam_Distr_Load+=['Dis_I='+str(Len_as)+'-(i-1)*'+str(Elem_Size)]

APDL_Beam_Distr_Load+=['Dis_J=(i)*'+str(Elem_Size)+'-'+str(Len_bs)]

APDL_Beam_Distr_Load+=["node_1=strcat('Node_',chrval(i))"]

APDL_Beam_Distr_Load+=["node_2=strcat('Node_',chrval(i+1))"]

APDL_Beam_Distr_Load+=['cmsel,s,node_1,node']

APDL_Beam_Distr_Load+=['cmsel,a,node_2,node']

APDL_Beam_Distr_Load+=['esln,s,1,corner']

APDL_Beam_Distr_Load+=['sfbeam,all,3,pres,'+str(Load_a_x)+','+str(Load_b_x)+',,,Dis_I,Dis_J']

殼單元均布荷載的施加方法與梁單元相似,重點在于單元坐標系z軸的確定。圖2為默認情況下,殼單元的單元坐標系。e1、e2、e3分別為單元坐標1軸、2軸、3軸;S1和S2為第一和第二參考方向;關鍵項KEYOPT(11)為單元x軸的默認方向,若KEYOPT(11)=0,則單元x軸的默認方向與質(zhì)心的第一個參數(shù)方向S1一致。

圖2 ANSYS中Shell181單元的默認單元坐標系Fig.2 Element coordinate system of Shell181 in ANSYS

參考方向S1和S2與單元的形函數(shù)有關。根據(jù)文獻[22],參考方向S1的計算方法如下

(1)

(2)

其中,{x}I、{x}J、{x}K、{x}L為總體坐標系下節(jié)點的坐標。

默認單元坐標軸e3由下式計算得到

e3=S1×S2/|S1×S2|

(3)

面劃分網(wǎng)格后,殼單元四個節(jié)點(IJKL)僅能在ANSYS程序中提取分析。因此,參考方向與單元坐標的確定通過APDL命令流實現(xiàn),其轉(zhuǎn)化代碼如下。首先,確定施加荷載的面(Area_Name),挑選其上編號最大的殼單元,定義四個節(jié)點的數(shù)組以分別存放節(jié)點坐標,此處以I節(jié)點為例。

APDL_Shell_Unif_Load+=['*dim,Area'+Area_Name+'_NodeI_Vec,array,3,1']

APDL_Shell_Unif_Load+=['Area'+Area_Name+'_NodeI_Vec(1,1)=nx(Area'+Area_Name+'_NodeI)']

APDL_Shell_Unif_Load+=['Area'+Area_Name+'_NodeI_Vec(2,1)=ny(Area'+Area_Name+'_NodeI)']

APDL_Shell_Unif_Load+=['Area'+Area_Name+'_NodeI_Vec(3,1)=nz(Area'+Area_Name+'_NodeI)']

然后,根據(jù)式(1)和式(2),計算得到殼單元的參考方向S1,S2與之相似。

APDL_Shell_Unif_Load+=['*voper,Area'+Area_Name+'_Parr1,Area'+Area_Name+'_NodeJ_Vec,sub,Area'+Area_Name+'_NodeI_Vec']

APDL_Shell_Unif_Load+=['*voper,Area'+Area_Name+'_Parr2,Area'+Area_Name+'_NodeK_Vec,sub,Area'+Area_Name+'_NodeL_Vec']

APDL_Shell_Unif_Load+=['*voper,Area'+Area_Name+'_Parr3,Area'+Area_Name+'_Parr1,add,Area'+Area_Name+'_Parr2']

APDL_Shell_Unif_Load+=['*voper,Area'+Area_Name+'_S1_A,Area'+Area_Name+'_Parr3,mult,1/4']

APDL_Shell_Unif_Load+=['Area'+Area_Name+'_S1_B=(Area'+Area_Name+'_S1_A(1,1)**2+Area'+Area_Name+'_S1_A(2,1)**2+Area'+Area_Name+'_S1_A(3,1)**2)**(1/2)']

APDL_Shell_Unif_Load+=['*voper,Area'+Area_Name+'_S1,Area'+Area_Name+'_S1_A,div,Area'+Area_Name+'_S1_B']

最后,根據(jù)式(3),計算得到殼單元坐標軸e3。得到殼單元的面外法向坐標軸,即可通過與梁單元相似的方法,編譯施加殼單元均布荷載的代碼。

APDL_Shell_Unif_Load+=['*voper,Area'+Area_Name+'_Parr7,Area'+Area_Name+'_S1,cross,Area'+Area_Name+'_S2']

APDL_Shell_Unif_Load+=['Area'+Area_Name+'_Parr8=(Area'+Area_Name+'_Parr7(1,1)**2+Area'+Area_Name+'_Parr7(1,2)**2+Area'+Area_Name+'_Parr7(1,3)**2)**(1/2)']

APDL_Shell_Unif_Load+=['*voper,Area'+Area_Name+'_Parr9,Area'+Area_Name+'_Parr9,div,Area'+Area_Name+'_Parr8']

APDL_Shell_Unif_Load+=['*dim,Area'+Area_Name+'_esys3,array,1,3']

APDL_Shell_Unif_Load+=['Area'+Area_Name+'_esys3(1,1)=Area'+Area_Name+'_Parr9(1,1)']

APDL_Shell_Unif_Load+=['Area'+Area_Name+'_esys3(1,2)=Area'+Area_Name+'_Parr9(2,1)']

APDL_Shell_Unif_Load+=['Area'+Area_Name+'_esys3(1,3)=Area'+Area_Name+'_Parr9(3,1)']

2 STAMT程序算例驗證

2.1 STAMT程序操作界面

STAMT程序的操作界面包括轉(zhuǎn)化程序(圖3(a))和使用說明(圖3(b))兩個基本模塊。

圖3 STAMT程序操作界面Fig.3 STAMT program operation interface

轉(zhuǎn)化程序模塊提供了梁單元劃分數(shù)的輸入、SAP2000文件錄入、APDL文件輸出、轉(zhuǎn)化進度顯示等功能。使用時,根據(jù)需要輸入梁單元網(wǎng)格劃分個數(shù),并點擊“瀏覽”按鈕選取轉(zhuǎn)化和導出的文件,或在輸入框中自定義工作目錄。完成上述步驟后,點擊“轉(zhuǎn)化”按鈕,就會出現(xiàn)轉(zhuǎn)化進度條,出現(xiàn)“轉(zhuǎn)化完成”時,代表模型轉(zhuǎn)化結束。

2.2 算例一:單層工業(yè)廠房

廠房基本布局為:縱向3跨,單跨4 m;橫向1跨,跨度12 m;屋架下弦高度6 m;檐口高度7.5 m;屋脊高度8.5 m。柱采用箱形鋼截面,梁采用工字形鋼截面,桁架采用工字形鋼桁架,屋面板為混凝土薄殼板。該結構的SAP2000模型如圖4(a)所示,轉(zhuǎn)換后的ANSYS模型如圖4(b)。

圖4 單層工業(yè)廠房模型Fig.4 Single-story industrial plant model

通過計算,對比兩個軟件的模態(tài)變形,如圖5、圖6所示。結果表明,SAP2000計算的前3階模態(tài)與ANSYS的計算結果基本一致。

圖5 未釋放自由度單層工業(yè)廠房模型的模態(tài)對比Fig.5 Modal comparison of single-storey industrial models with unreleased degree of freedom

圖6 釋放自由度后單層工業(yè)廠房模型的模態(tài)對比Fig.6 Modal comparison of single-storey industrial models with released degree of freedom

為更加準確地驗證STAMT程序的轉(zhuǎn)化效果,分別按釋放桁架自由度和不釋放桁架自由度兩種情況,對比SAP2000模型和ANSYS模型的質(zhì)量和周期,如表2所示。可以看出,SAP2000模型的總質(zhì)量和ANSYS模型一致。對于不釋放自由度的模型,SAP2000的計算結果與ANSYS程序相近,誤差平均值為2.44%,最大誤差不超過5%。對于釋放自由度后的模型,SAP2000計算的前3階周期與ANSYS相近。

表2 單層工業(yè)廠房周期對比Table 2 Period comparison of single-storey industrial plant

但在較大周期(>4階)時,二者的計算誤差增加,這是由于SAP2000軟件使用的是節(jié)點集中質(zhì)量,而ANSYS軟件在釋放梁端自由度之后,模態(tài)計算無法打開集中質(zhì)量開關(lumpm, on),從而造成二者的計算結果出現(xiàn)差異,且誤差在高階模態(tài)下不斷增加。

整體而言,通過STAMT程序可以實現(xiàn)SAP2000模型向ANSYS模型的準確轉(zhuǎn)化,特別是低階模態(tài)的計算誤差不超過5%。然而,對于釋放自由度后,集中質(zhì)量與一致質(zhì)量的轉(zhuǎn)化問題,仍需在今后的軟件優(yōu)化中進一步研究。

2.3 算例二:多層框架結構

框架結構長30 m,寬18 m,高18 m,結構由下部鋼筋混凝土構件組成。梁柱構件截面種類為矩形或方形。整體框架結構的SAP2000模型如圖7(a)所示。通過STAMT程序轉(zhuǎn)化SAP2000模型,生成APDL命令流,將命令流導入ANSYS中,創(chuàng)建的ANSYS有限元模型如圖7(b)所示。

圖7 多層框架結構模型Fig.7 Multi-storey framework structural model

對比兩個模型的模態(tài)變形,如圖8所示?？梢钥闯?SAP2000計算的模態(tài)變形圖與ANSYS的計算結果一致,表明轉(zhuǎn)化后的模型滿足要求。

圖8 多層框架結構模型的模態(tài)對比Fig.8 Modal comparison of multi-storey frame structures

進一步比較二者的質(zhì)量以及計算的前5階周期,如表3所示。結果表明,兩個模型的結構總質(zhì)量一致,而周期的誤差在2%以內(nèi),且平均誤差為1.32%,說明轉(zhuǎn)化模型的計算精度滿足要求。

表3 多層框架結構質(zhì)量及周期對比Table 3 Comparison of weight and period for multi-storey frame structures

通過以上的算例分析,證明了本文所編譯開發(fā)的STAMT程序可以快速、準確地將SAP2000模型轉(zhuǎn)化為ANSYS有限元分析模型,可滿足相關技術人員對模型轉(zhuǎn)化的需求。

3 結論與展望

基于Python語言,編譯并開發(fā)了STAMT程序,實現(xiàn)了SAP2000結構設計分析模型向ANSYS有限元模型的轉(zhuǎn)化,并通過兩個算例的計算驗證了該程序的可行性及準確性。相比于現(xiàn)有的模型轉(zhuǎn)化程序,本文開發(fā)的STAMT程序具有以下優(yōu)點。

(1)豐富了單元類型、材料屬性、截面形式、荷載形式的轉(zhuǎn)化類型,提高了模型轉(zhuǎn)化的精確性。

(2)涵蓋了節(jié)點坐標系變換、梁單元坐標系變換、梁端自由度釋放、創(chuàng)建剛性域等必備功能,進一步滿足了模型轉(zhuǎn)化過程中的功能性需求。

(3)以ANSYS的關鍵點(Keypoint)對應SAP2000的結點(Joint),使用先建立幾何模型后生成有限元模型的思路,實現(xiàn)了梁單元網(wǎng)格數(shù)量的自定義設置以及板殼荷載向梁的傳遞。

開發(fā)的STAMT程序可為相關從業(yè)人員更加快速、準確地實現(xiàn)模型信息的轉(zhuǎn)化,為結構分析結果的對比驗證提供便捷。但由于結構設計分析軟件與通用有限元軟件之間的差異性,難以實現(xiàn)SAP2000模型向ANSYS模型的完全轉(zhuǎn)化,STAMT程序中存在的不足會隨著各分析軟件的不斷更新以及開發(fā)者的修復改進,而在將來的研究中予以完善。