一個(gè)金屬長(zhǎng)方體，我們需要對(duì)其做拉伸的加載約束示意圖如圖1，并在完成后采用Python命令流讀取參考點(diǎn)的位移、體積、應(yīng)變隨加載時(shí)間的變化情況。

                                                 圖1 金屬長(zhǎng)方體約束加載示意圖 

對(duì)于圖上的命令流，我們分為幾部分拆開研究，首先是下面第一部分需要在Abaqus導(dǎo)入相關(guān)Python庫(kù)的代碼，比如matplotlib庫(kù)、numpy庫(kù)、xlwt庫(kù)和math庫(kù)。

#! /user/bin/python
# -*- coding:UTF-8 -*-
import matplotlib.pyplot as plt
import numpy as np
from visualization import *
from odbAccess import *
import xlwt
import math

• 第3行導(dǎo)入Python經(jīng)常用于繪圖的matplotlib庫(kù)；
  

• 第4行導(dǎo)入Python處理數(shù)組矩陣的numpy庫(kù)；

• 第5行和第6行用于導(dǎo)入Abaqus的后處理功能，即visualization和odbAccess；

• 第7行用于導(dǎo)入Python與excel進(jìn)行交互的xlwt庫(kù)；

• 第8行用于導(dǎo)入Python的math庫(kù)。

在利用Python生成excel數(shù)據(jù)之前，先采用一個(gè)getInputs函數(shù)（代碼如下）生成圖2所示的對(duì)話框與用戶進(jìn)行交互，需要用戶在在其中輸入模型的名字（model name），部件例子的名字（instance name）以及后處理odb的名字（odbname）。這樣做有個(gè)好處，每次都可以根據(jù)不同模型、不同部件實(shí)例和相應(yīng)的后處理odb名字進(jìn)行相應(yīng)結(jié)果提取。

modelName,instanceName,odbname= getInputs( 
 fields=(('Model Name:', 'Test'),('Instance Name:', 'Part-2-1'),('odbname:', 'TEST.odb')), 
  label='Enter information',  dialogTitle='Enter information.')

                                                    圖2 getInputs函數(shù)生成的對(duì)話框

在計(jì)算體積和各邊應(yīng)變之前，需要計(jì)算長(zhǎng)方體不同邊長(zhǎng)即長(zhǎng)、寬、高的初始長(zhǎng)度。代碼如下：

node = mdb.models[modelName].rootAssembly.instances[instanceName].nodes
Xmin = 9999
Xmax = -9999
Ymin = 9999
Ymax = -9999
Zmin = 9999
Zmax = -9999
for i in range(len(node)):   
    x = node[i].coordinates[0]    
    y = node[i].coordinates[1]    
    z = node[i].coordinates[2]    
    if Xmin > x:        
        Xmin = x    
    elif Xmax < x:        
        Xmax = x            
    if Ymin > y:        
        Ymin = y    
    elif Ymax < y:        
        Ymax = y            
    if Zmin > z:        
        Zmin = z        
        continue    
    if Zmax < z:        
        Zmax = z        
        continue
print 'Xmin,Xmax,Ymin,Ymax,Zmin,Zmax=',(Xmin, Xmax, Ymin, Ymax, Zmin, Zmax)
CubeWidth=abs(Xmax-Xmin)  
# the length of model along x directionCubeHeight=abs(Ymax-Ymin) 
# the length of model along y directionCubeLength=abs(Zmax-Zmin) 
# the length of model along z direction

• 第1行基于之前輸入的modelname，instancename來(lái)調(diào)用部件實(shí)例（instance）的節(jié)點(diǎn)。

• 第8行for i in range(len(node))：表示對(duì)每個(gè)節(jié)點(diǎn)進(jìn)行循環(huán)，每次循環(huán)中都要提取該節(jié)點(diǎn)的x，y，z坐標(biāo)信息。
  

• Xmin，Xmax為通過(guò)循環(huán)找出X方向的最小、最大X坐標(biāo)，Ymin、Ymax以及Zmin、Zmax都同理。計(jì)算這部分的代碼可能大家表面無(wú)法直接看懂理解，我舉一個(gè)例子好了。

例子：對(duì)于第12行到第15行，比如當(dāng)i=0時(shí)，隨便提取了第1個(gè)節(jié)點(diǎn)的x，y，z信息（0,0.1,0.2），由于Xmin=9999>x=0，所以x=0賦予Xmin；同樣的，當(dāng)i=1時(shí)，又隨便提取模型的一個(gè)節(jié)點(diǎn)(0.05,0.1,0.2)，由于現(xiàn)在Xmin=0<x=0.05，則執(zhí)行Xmax=-9999<x=0.05，所以x=0.05賦予Xmax，而Xmin不變；當(dāng)i=2時(shí)，隨便提取一個(gè)節(jié)點(diǎn)(0.07,0.1,0.2)，注意到此時(shí)，Xmin=0<x=0.07并且Xmax=0.05<x=0.07，因此x=0.07賦予Xmax，即此時(shí)的最大值賦予給Xmax，而Xmin不變；而如果i循環(huán)到節(jié)點(diǎn) (0.1,0.1,0.2)，注意到此時(shí)，Xmin=0<x=0.1并且Xmax=0.05<x=0.1，因此x=0.1賦予Xmax，即此時(shí)的最大值賦予給Xmax，而Xmin不變；由于x=0是模型x方向的最小值，x=0.1是模型x方向的最大值，因此如果當(dāng)i循環(huán)到又一個(gè)節(jié)點(diǎn)(0.06,0.1,0.2)，此時(shí)，Xmin=0<x=0.06并且Xmax=0.1>x=0.06，則Xmin和Xmax則處于不變，因?yàn)閤在這兩者之間，因此通過(guò)以上程序就可以找出模型在x方向的最大值和最小值（結(jié)果見圖3），然后最大值減去最小值則得到了模型沿x方向的長(zhǎng)度。同樣的方式得到模型沿y方向和z方向的長(zhǎng)度。

圖3 查看Xmin、Xmax、Ymin、Ymax、Zmin、Zmax

另一方面，我們可以如圖4進(jìn)行相應(yīng)驗(yàn)證，我們先查看下模型在x方向的寬度，在Abaqus中的query中選取distance，然后選中模型中的兩個(gè)點(diǎn)，可以看到這兩個(gè)點(diǎn)的坐標(biāo)是point1（0,0.1,0.2），point2（0.1,0.1,0.2），則兩點(diǎn)之間在x方向的距離為0.1，即為模型的寬度（模型沿x方向的長(zhǎng)度）。Abaqus模型計(jì)算的結(jié)果與之前程序計(jì)算的結(jié)果相互驗(yàn)證了。

圖4 查看模型在x方向的寬度

odb=openOdb(path=odbname)
wbkName='NodalDisplacement'
wbk1=xlwt.Workbook()
sheet=wbk1.add_sheet('Sheet1')
RefPointSet=odb.rootAssembly.nodeSets['U-XYZ']
frameRepository=odb.steps['Step-1'].frames
for i in range(3):    
    sheet.write(0,i,('U'+str(i+1)+'(U-XYZ)'))

• 第1行表示打開odb計(jì)算結(jié)果文件。

• 第2行表示生成的excel文件名為NodalDisplacement。

• 第3行和第4行表示添加Sheet1工作表。

• 第5行建立節(jié)點(diǎn)SetU-XYZ，方便后續(xù)對(duì)該Set進(jìn)行位移提取。

• 第6行表示Step-1里的幀(frames)對(duì)象。

• 第7行和第8行通過(guò)循環(huán)在excel第1行第i列（也就是代碼中的0，i），寫入U(xiǎn)1(U-XYZ)，U2(U-XYZ)，U3(U-XYZ)，如圖5所示，其中Python編程用+可以實(shí)現(xiàn)字符連綴，str(i)可以將數(shù)字變量i改為字符i。

iframes=len(frameRepository)
RefDataU1= np.zeros((iframes,1))
RefDataU2= np.zeros((iframes,1))
RefDataU3= np.zeros((iframes,1))
jingdu=15
for i in range(iframes):    
    U=frameRepository[i].fieldOutputs['U']    
    RefU=U.getSubset(region=RefPointSet)   
    RefUValues=RefU.values    
    RefDataU1[i]=RefUValues[0].data[0]   
    sheet.write(i+1,0,round(RefDataU1[i],jingdu))   
    RefDataU2[i]=RefUValues[0].data[1]    
    sheet.write(i+1,1,round(RefDataU2[i],jingdu)) 
    RefDataU3[i]=RefUValues[0].data[2]   
    sheet.write(i+1,2,round(RefDataU3[i],jingdu))
wbk1.save(wbkName+'.xls')

• 第一行是統(tǒng)計(jì)frames的幀數(shù)。

• 第2行到第4行基于numpy定義三個(gè)數(shù)組存儲(chǔ)三個(gè)方向的位移，并事先賦予初值零。

• 第5行設(shè)置一個(gè)變量，方便用于控制輸出到excel里數(shù)的精度。

• 第6行到第15行就是通過(guò)for循環(huán)，一幀一幀地輸出U-XYZ這個(gè)點(diǎn)Set的三個(gè)方向的位移值，并存儲(chǔ)到excel。U-XYZ這個(gè)點(diǎn)Set如圖5所示。

圖5 U-XYZ點(diǎn)Set

• 對(duì)于第7行表示提取該幀的位移U的場(chǎng)變量并賦值給U，然后在第8行采用getSubset基于U提取RefPointSet區(qū)域的位移值。

• 圖6在ABAQUS命令行接口（基于Python的Abaqus靜力分析操作實(shí)例）表示的應(yīng)該是最后一次循環(huán)后各變量的值，可以看出RefU表示的nodeSets['U-XYZ']的位移值，圖6中顯示為一個(gè)索引；而RefUValues則表示該節(jié)點(diǎn)集中每個(gè)節(jié)點(diǎn)的位移值，體現(xiàn)在values中；而values[0]表示第一個(gè)節(jié)點(diǎn)的位移值，由于存在三個(gè)方向的位移，U1,U2,U3,則后續(xù)通過(guò)data[0]，data[1]，以及data[2]分別表示三個(gè)方向的位移值。
  

圖6 運(yùn)行相關(guān)結(jié)果查看（通過(guò)ABAQUS命令行接口）

特別的，由于U-XYZ只有一個(gè)點(diǎn)存在于Set中，如果我們?cè)谶@里用values[1]，看會(huì)出現(xiàn)什么樣的效果？如圖7。可以明顯看到Sequence index out of range表示超出了范圍，因?yàn)镽efUValues[1]表示提取第二個(gè)節(jié)點(diǎn)，但我們這里只有一個(gè)節(jié)點(diǎn)存在于U-XYZ點(diǎn)Set中，所以出現(xiàn)了范圍超過(guò)的報(bào)錯(cuò)提示。

圖7 可以存在RefUValues[1]嗎？（通過(guò)ABAQUS命令行接口）

需要說(shuō)一下的，如果大家覺得ABAQUS命令行接口的查看范圍太窄，可以往上拉動(dòng)窗口，這樣就可以看到更多的代碼，如圖8所示，但與之相反的就是模型就會(huì)顯得很小，與我們一般的常規(guī)abaqus操作有點(diǎn)相反。

圖8 ABAQUS命令行接口顯示更多的內(nèi)容

再對(duì)代碼中的Python round函數(shù)做個(gè)說(shuō)明，通過(guò)下面的例子我們可以看到round函數(shù)是四舍五入函數(shù)，一般為四舍五入后為整數(shù)，round(a,b)中的a為需要進(jìn)行四舍五入的數(shù)，而b為需要保留的小數(shù)位數(shù)，如圖9所示。

圖9 關(guān)于Python round 函數(shù)的例子

最后就是第17行代碼，這個(gè)代碼也很重要，不輸入這個(gè)代碼就不會(huì)生成excel表，一般生成excel位置為設(shè)置好的工作目錄（比如D:\temp），生成好的excel如圖10所示。而excel總共有625行，說(shuō)明一共有625幀，也就是上述代碼的變量iframes為625。

圖10 生成的NodalDisplacement.xls的excel文件

04

結(jié)果驗(yàn)證

將上述excel表格中的結(jié)果繪制成曲線，即U-XYZ點(diǎn)三個(gè)方向的位移隨幀數(shù)變化的曲線，如圖11所示。根據(jù)圖12，再結(jié)合odb文件中的變形圖（第623幀）可以判斷U-XYZ點(diǎn)在后期的X(U1)為負(fù)，Y(U2)為正，Z(U3)為正，與圖11曲線是一致的。

圖11 U-XYZ點(diǎn)三個(gè)方向的位移隨幀數(shù)變化的曲線

圖12 模型后處理結(jié)果圖（第623幀）

最后再與在Abaqus里進(jìn)行GUI（即點(diǎn)點(diǎn)點(diǎn)）對(duì)U-XYZ三個(gè)方向位移處理的結(jié)果圖進(jìn)行對(duì)比，即通過(guò)XY Data對(duì)U-XYZ三個(gè)方向的位移進(jìn)行輸出繪圖，注意到坐標(biāo)軸為時(shí)間，與Python后處理圖11的幀數(shù)是一一對(duì)應(yīng)的。可以知道，Python后處理結(jié)果與ABAQUS的GUI操作得到的結(jié)果幾乎完全一致。

圖13 Abaqus里進(jìn)行GUI操作的后處理結(jié)果（關(guān)于U-XYZ三方向位移）

該模型關(guān)于體積、應(yīng)變、應(yīng)力等方面（包括場(chǎng)輸出和歷史輸出）的Python導(dǎo)入excel的后處理留在下期，敬請(qǐng)關(guān)注。

最后，歡迎大家通過(guò)微信公眾號(hào)聯(lián)系我們。

微信公眾號(hào)：320科技工作室。

同時(shí)，更多ABAQUS二次開發(fā)的精彩內(nèi)容，歡迎大家關(guān)注微信公眾號(hào)： 土木科研編程與數(shù)值模擬