//Logo Image
作者: 許銘修(2002-08-22);核可:徐業良(2002-08-22)

建立程式化之ANSYS輸入檔

在結構最佳化設計過程中,常需整合最佳化設計程式與有限元素分析軟體,反覆執行有限元素分析,以提供最佳化設計成是所需之數值資訊。為了使整個設計迭代過程為一自動化過程,因此如何將最佳化設計程式以及有限分析軟體所需之輸入檔予以程式化,並且讀取分析計算後之結果檔,對於此類最佳化設計程序而言是相當重要的技術。本文即在介紹如何利用MATLAB將一ANSYS之輸入檔予以程式化成一副程式,使其可以以批次(batch)方式進行有限元素分析。

1. ANSYS輸入檔

結構形態最佳化設計(shape optimization)的目的是希望能找到一組最佳的結構邊界形態,使得結構的剛性為最大而又不會產生破壞。一般的作法是以spline曲線來表示要設計的結構邊界,而設計變數則為此spline曲線的控制點的法線方向移動量。圖1所示為一結構形態最佳化設計範例之幾何形狀與邊界條件,其中Keypoint 111spline曲線的控制點,其餘各點則以直線連接。於ANSYS有限元素分析軟體中,可以使用GUI介面進行模型的建立與分析,也可以選擇將所要進行的工作指令寫成輸入檔,以批次方式進行分析。為了使最佳化設計工作能自動迭代,我們選擇撰寫分析輸入檔,表1所示為圖1中之形態最佳化設計範例之ANSYS輸入檔。

1. 形態最佳化設計範例之幾何形狀與邊界條件

1. ANSYS輸入檔

!INITIAL DEFINITION

        /TITLE,2D SHAPE OPTIMIZATION EXAMPLE

        /PREP7

        ANTYPE,STATIC

 

        ET,1,PLANE82

        KEYOPT,1,3,3

        R,1,1

        MP,EX,1,210000

        MP,NUXY,1,0.3

        MP,DENS,1,1

 

!INPUT KEYPOINTS

        K,1,1560.000000,320.000000,0.000000

        K,2,1540.000000,300.000000,0.000000

        K,3,1500.000000,300.000000,0.000000

        K,4,1439.992000,279.996700,0.000000

        K,5,1160.041000,280.005000,0.000000

        K,6,920.000500,279.997500,0.000000

        K,7,699.989700,280.015100,0.000000

        K,8,600.006800,300.005700,0.000000

        K,9,460.003500,299.998500,0.000000

        K,10,400.000400,319.999800,0.000000

        K,11,360.000000,360.000000,0.000000

        K,12,0.000000,0.000000,0.000000

        K,13,360.000000,720.000000,0.000000

        K,14,360.000000,0.000000,0.000000

        K,15,0.000000,720.000000,0.000000

        K,16,1560.000000,0.000000,0.000000

        K,17,1920.000000,0.000000,0.000000

        K,18,1920.000000,360.000000,0.000000

        K,19,1560.000000,360.000000,0.000000

 

!GENERATING LINES

        FLST,2,11,3

        FITEM,2,1  

        FITEM,2,2  

        FITEM,2,3  

        FITEM,2,4  

        FITEM,2,5  

        FITEM,2,6  

        FITEM,2,7  

        FITEM,2,8  

        FITEM,2,9  

        FITEM,2,10 

        FITEM,2,11 

        SPLINE,P51X

        LSTR,      11,      13 

        LSTR,      13,      15 

        LSTR,      15,      12 

        LSTR,      12,      14 

        LSTR,      14,      16 

        LSTR,      16,      17 

        LSTR,      17,      18 

        LSTR,      18,      19 

        LSTR,      19,       1 

 

!GENERATING AREA

        AL,ALL

 

!MESH

        ESIZE,15

        AMESH,ALL

 

!APPLY BOUNDARY CONDITIONS

        NSEL,S,LOC,X,0

        F,ALL,FX,-5000/97

        NSEL,ALL

        NSEL,S,LOC,X,1920,1920

        F,ALL,FX,5000/49

        NSEL,ALL

        NSEL,S,LOC,Y,0

        D,ALL,UY

        NSEL,ALL

        /OUTPUT,AREA,OUT,D:\research\shape_2d\Example_01\

        ASUM

        /OUTPUT,TERM

        FINISH

 

!SOLVE

        /SOLU

        EQSLV,PCG

        SOLVE

        FINISH

 

!POSTPROCESS

        /POST1

        /OUTPUT,STRESS,OUT,D:\research\shape_2d\Example_01\

        NSEL,S,F,F,-100000,100000

        PRNSOL,U

        NSEL,ALL

        PRNSOL,S,PRIN

        /OUTPUT,TERM

        FINISH

        /CLEAR,NOSTART

        /EXIT

2 設計變數的處理

結構形態最佳化設計主要在找尋一組最佳的邊界形態,就2-D結構而言,即是找尋一組最佳的邊界曲線,而最常使用的曲線為spline曲線。Spline曲線是由一組控制點所決定而成,移動任一控制點的位置,即會獲得不同的曲線形態。因此結構形態最佳化設計的問題,即可轉變為在給定的初始控制點位置下,找尋一組最佳的控制點移動量,亦即我們可以設定設計變數為每一控制點於法線方向之移動量。

在初始設計的控制點位置以及法線方向角度為已知的情況下,新的控制點位置可由式(1)與式(2)獲得:

                                                                                                (1)

                                                                                                (2)

其中xiyi表示第i個控制點移動後的xy座標,x0y0表示初始控制點的xy座標,di為第i個控制點於法線方向之移動量,而qii個控制點的法線方向與x軸的夾角角度。因此程式化後之ANSYS輸入檔程式的要求是只需輸入di即可自動產生完整之ANSYS輸入檔,以下及逐步介紹建立程式化後之ANSYS輸入檔所需指令。

3. ANSYS輸入檔程式化所需指令介紹

function 指定程式為函式型程式

ANSYS輸入檔程式指定為函式(function)型程式,我們即可將數字資料矩陣輕易的輸入程式中,因為函式型程式可以接受輸入變數,並且也可以將計算結果輸出。指定程式為函式型程式的方法是在程式起始處使用“function”指令,如表2所示即為一形態最佳化設計程式之範例,function指令後跟著為該程式的檔案名稱,此處之檔案名稱為Example_inputfile,接著給予函式變數名稱,此處共有兩個變數輸入,分別為產生之檔案路徑(work_path),以及控制點法線方向移動量(delta_input),輸入函式之變數只須要是矩陣形式即可。控制點法線方向移動量輸入程式後,即可依式(1)與式(2)計算出新的控制點位置。

2. 指定程式為函式型程式

function Example_inputfile(work_path,delta_input)

 

fopen 開啟檔案指令

fopen指令可用於開啟已存在之檔案,或是建立一新檔案,表3所示為一使用範例,fem_inp為一指標變數,用以表示該開啟檔案,指標變數名稱可以字以給定。Fopen中包含兩個輸入變數,首先必須給予要開啟檔案的路徑、檔名與副檔名,如表3中所示表示將開啟D:\shape_2d\fem_inp.inp檔案。接著給予工作模式,亦即指定是要讀取檔案(r),或是撰寫檔案(w),此處是要產生一ANSYS輸入檔,因此選擇使用“w”,又由於所要寫的檔案為文字檔,因此加一“t”

3. fopen指令使用範例

fem_inp=fopen('D:\shape_2d\fem_inp.inp','wt');

 

fprintf 將指定的文字、數字印於檔案上

檔案開啟後即可在檔案中寫上東西,fprintf指令即可將指定的字串或浮點印於檔案上,表4所示為fprintf指令之使用範例。fprintf指令變數有三部分,首先是指定要將字串印上的檔案指標,此指標即為在fopen指令中所指定的檔案指標,如表4中的範例,指定要印上的檔案為fem_inp指標所表示的檔案。

第二個變數部分為指定要使用的格式,常用的格式有“%s”表示為字串,“%f”表示為浮點,“\n”表示換行。第二數部分的給定與第三變數部分有關,第三部分為所要印上的字串或浮點。以表4中的範例來說,fprintf(fem_inp,'%s\n','!INITIAL DEFINITION')的第三變數部分為一字串'!INITIAL DEFINITION',因此第二變數部分中就必須使用“%s”,因為希望在印完字串後換行,因此多加了“\n”

較複雜的範例如表4的最後一行指令所示,此指令中所要印出的內容包含有字串與浮點兩種,其中'K,'以及兩個','是字串的部分,變數keypoint_matrix(kpoint_i,1)x_coor以及y_coor為浮點的部分,由於相互間隔,所以第二變數部分中的給定將是'%s%f%s%f%s%f'方式,又由於希望keypoint_matrix(kpoint_i,1)所表示的浮點數不要顯示出小數點後的位數,因此可以在%f中加入.0”成為“%.0f,加上換行的動作,所以整個格式的給定成為'%s%.0f%s%f%s%f\n'

於第三部分所要印出的字串或浮點之間,以“,”來作區隔,只要格式的給定正確,第三部分的內容可以一直增加,然而除非必要,不然不建議將所有要印出的文字或數字在一個指令中完成,因為這樣會增加閱讀的複雜度,使得程式撰寫時容易出錯。假設keypoint_matrix(kpoint_i,1)=1x_coor=10y_coor=10,則表4中最後一行指令的結果將為:

        K,1,10.000000,10.000000

4. fprintf指令使用範例

fprintf(fem_inp,'%s\n','!INITIAL DEFINITION');    

fprintf(fem_inp,'%s\n','/TITLE, 2D SHAPE OPTIMIZATION EXAMPLE ');

fprintf(fem_inp,'%s\n','      /PREP7');

 

fprintf(fem_inp,'%s%.0f%s%f%s%f\n','K,',keypoint_matrix(kpoint_i,1),',',x_coor,',',y_coor);

 

fclose 關閉檔案指令

當完成檔案時可以使用fclose指令將檔案關閉,fclose指令的變數即是要關閉的檔案的指標,如表5範例中所示,要關閉的檔案是稱為fem_inp指標的檔案,亦即於表3中所開啟的檔案。

有時因為須要而在程式中同時開啟數個檔案,此時若是一個一個關閉就太過麻煩了,因此可以將指令變數中的檔案指標改為如表5第二行指令所示的'all',即可將所有檔案一次關閉。

5. fclose指令使用範例

fclose(fem_inp)

fclose('all')

 

4. 程式範例

6所示為由表1中之ANSYS輸入檔所改寫而來的副程式,此副程式之輸入變數有work_path,表示工作路徑,所產生之ANSYS輸入檔,以及ANSYS分析後的輸出檔都將存放於此路徑下。另一輸入變數為delta_input,此變數即為各個控制點於法線方線的移動量。若是令work_path='D:\research\shape_2d\Example_01\',而delta_input=[0,0,0,0,0,0,0,0,0,0,0],則所得之ANSYS輸入檔就是表1中所示的輸入檔。

6. 程式範例

function Example_inputfile(work_path,delta_input)

%

%     Shape optimization example

%

%

%                     Ming-Hsiu Hsu

%                     Aug.,12,2002

 

        if work_path(length(work_path))~='\'

                work_path=[work_path,'\'];

        end;

 

 

        keypoint_matrix=[

                        1,1560.000,320.0000,0.000000

                        2,1540.000,300.0000,0.000000 

                        3,1500.000,300.0000,0.000000

                        4,1439.992,279.9967,0.000000

                        5,1160.041,280.0050,0.000000

                        6,920.0005,279.9975,0.000000 

                        7,699.9897,280.0151,0.000000

                        8,600.0068,300.0057,0.000000

                        9,460.0035,299.9985,0.000000

                        10,400.0004,319.9998,0.000000

                        11,360.0000,360.0000,0.000000

                        12,0.000000,0.000000,0.000000

                        13,360.0000,720.0000,0.000000

                        14,360.0000,0.000000,0.000000

                        15,0.000000,720.0000,0.000000

                        16,1560.000,0.000000,0.000000

                        17,1920.000,0.000000,0.000000

                        18,1920.000,360.0000,0.000000

                        19,1560.000,360.0000,0.000000

                        ];

 

        angle_matrix=[

                        1,180

                        2,67.5

                        3,80.7823

                        4,80.7831

                        5,90

                        6,90

                        7,84.3444

                        8,84.3482

                        9,80.7839

                        10,58.2825

                        11,0

                        ];

 

 

 

        eval(['fem_inp=fopen(''',work_path,'fem_inp.inp'',''wt'');']);

 

 

        fprintf(fem_inp,'%s\n','!INITIAL DEFINITION');    

        fprintf(fem_inp,'%s\n','/TITLE, 2D SHAPE OPTIMIZATION EXAMPLE ');

        fprintf(fem_inp,'%s\n','      /PREP7');

        fprintf(fem_inp,'%s\n','      ANTYPE,STATIC');

        fprintf(fem_inp,'\n');

        fprintf(fem_inp,'%s\n','      ET,1,PLANE82');

        fprintf(fem_inp,'%s\n','      KEYOPT,1,3,3');

        fprintf(fem_inp,'%s\n','      R,1,1');

        fprintf(fem_inp,'%s\n','      MP,EX,1,210000');

        fprintf(fem_inp,'%s\n','      MP,NUXY,1,0.3');

        fprintf(fem_inp,'%s\n','      MP,DENS,1,1');

        fprintf(fem_inp,'\n');

 

        fprintf(fem_inp,'%s\n','!INPUT KEYPOINT');

       

        [kpoint_n,kpoint_m]=size(keypoint_matrix);

        [angle_n,angle_m]=size(angle_matrix)

 

        for kpoint_i=1:kpoint_n

                if kpoint_i<=angle_n

                x_coor=keypoint_matrix(kpoint_i,2)

                        +delta_input(kpoint_i)*cos(angle_matrix(kpoint_i,2)*pi/180);

                y_coor=keypoint_matrix(kpoint_i,3)

                        +delta_input(kpoint_i)*sin(angle_matrix(kpoint_i,2)*pi/180);

                z_coor=keypoint_matrix(kpoint_i,4);

                else

                x_coor=keypoint_matrix(kpoint_i,2);

                y_coor=keypoint_matrix(kpoint_i,3);

                z_coor=keypoint_matrix(kpoint_i,4);                       

                end;

 

                fprintf(fem_inp,'%s%.0f%s%f%s%f%s%f\n','    K,', …

                        keypoint_matrix(kpoint_i,1),',',x_coor,',',y_coor,',',z_coor);

        end;

 

        fprintf(fem_inp,'\n%s\n','!CREATING LINES');

        fprintf(fem_inp,'%s\n','      FLST,2,11,3 ');

        fprintf(fem_inp,'%s\n','      FITEM,2,1   ');

        fprintf(fem_inp,'%s\n','      FITEM,2,2   ');

        fprintf(fem_inp,'%s\n','      FITEM,2,3   ');

        fprintf(fem_inp,'%s\n','      FITEM,2,4   ');

        fprintf(fem_inp,'%s\n','      FITEM,2,5   ');

        fprintf(fem_inp,'%s\n','      FITEM,2,6   ');

        fprintf(fem_inp,'%s\n','      FITEM,2,7   ');

        fprintf(fem_inp,'%s\n','      FITEM,2,8   ');

        fprintf(fem_inp,'%s\n','      FITEM,2,9   ');

        fprintf(fem_inp,'%s\n','      FITEM,2,10  ');

        fprintf(fem_inp,'%s\n','      FITEM,2,11  ');

        fprintf(fem_inp,'%s\n','      SPLINE,P51X ');

        fprintf(fem_inp,'%s\n','      LSTR,      11,      13  ');

        fprintf(fem_inp,'%s\n','      LSTR,      13,      15  ');

        fprintf(fem_inp,'%s\n','      LSTR,      15,      12  ');

        fprintf(fem_inp,'%s\n','      LSTR,      12,      14  ');

        fprintf(fem_inp,'%s\n','      LSTR,      14,      16  ');

        fprintf(fem_inp,'%s\n','      LSTR,      16,      17  ');

        fprintf(fem_inp,'%s\n','      LSTR,      17,      18  ');

        fprintf(fem_inp,'%s\n','      LSTR,      18,      19  ');

        fprintf(fem_inp,'%s\n','      LSTR,      19,       1  ');

 

        fprintf(fem_inp,'\n%s\n','!CREATING AREA');

        fprintf(fem_inp,'%s\n','      AL,ALL');

 

        fprintf(fem_inp,'\n%s\n','!MESH');

        fprintf(fem_inp,'%s\n','      ESIZE,15');

        fprintf(fem_inp,'%s\n','      AMESH,ALL');

 

        fprintf(fem_inp,'\n%s\n','!APPLY BOUNDARY CONDITIONS');

        fprintf(fem_inp,'%s\n','      NSEL,S,LOC,X,0');

        fprintf(fem_inp,'%s\n','      F,ALL,FX,-5000/97');

        fprintf(fem_inp,'%s\n','      NSEL,ALL');

        fprintf(fem_inp,'%s\n','      NSEL,S,LOC,X,1920,1920');

        fprintf(fem_inp,'%s\n','      F,ALL,FX,5000/49');

        fprintf(fem_inp,'%s\n','      NSEL,ALL');

        fprintf(fem_inp,'%s\n','      NSEL,S,LOC,Y,0');

        fprintf(fem_inp,'%s\n','      D,ALL,UY');

        fprintf(fem_inp,'%s\n','      NSEL,ALL');

        eval(['fprintf(fem_inp,''%s\n'',''/OUTPUT,AREA,OUT,',work_path,''');']);

        fprintf(fem_inp,'%s\n','      ASUM');

        fprintf(fem_inp,'%s\n','      /OUTPUT,TERM');

        fprintf(fem_inp,'%s\n','      FINISH');

 

        fprintf(fem_inp,'\n%s\n','!SOLVE');

        fprintf(fem_inp,'%s\n','      /SOLU');

        fprintf(fem_inp,'%s\n','      EQSLV,PCG');

        fprintf(fem_inp,'%s\n','      SOLVE');

        fprintf(fem_inp,'%s\n','      FINISH');

 

        fprintf(fem_inp,'\n%s\n','!POSTPROCESS');

        fprintf(fem_inp,'%s\n','      /POST1');

        eval(['fprintf(fem_inp,''%s\n'',''/OUTPUT,STRESS,OUT,',work_path,''');']);

        fprintf(fem_inp,'%s\n','      NSEL,S,F,F,-100000,100000');

        fprintf(fem_inp,'%s\n','      PRNSOL,U');

        fprintf(fem_inp,'%s\n','      NSEL,ALL');

        fprintf(fem_inp,'%s\n','      PRNSOL,S,PRIN');

        fprintf(fem_inp,'%s\n','      /OUTPUT,TERM'); 

        fprintf(fem_inp,'%s\n','      FINISH');

        fprintf(fem_inp,'%s\n','      /CLEAR,NOSTART');

        fprintf(fem_inp,'%s\n','      /EXIT');

 

        fclose('all');