!************************* ***橋有限元模型 ******************************

!**橋的整體有限元模型建立。該模型為魚脊骨模型。先按照梁的自然節(jié)點（拉索位置）
!和主塔的自然節(jié)點（截面變化處、縱橫梁交叉點）建立三維梁單元，然后分別讀入前面
!建立的梁、塔的截面信息，生成空間整體三維有限元梁柱子體系。為了后面建立索的方便，
!索的兩段的節(jié)點必須是按照等差數(shù)列的方式排列的。并且用參數(shù)獲得索的第一個節(jié)點號、
!最后一個節(jié)點號，以及節(jié)點增加步長。然后使用循環(huán)語句分別生成八排索。

!**************************** 本程序的單位制 *******************************

!質(zhì)量；噸； 力：牛頓； 長度：毫米； 時間：秒； 應力及彈性模量：兆帕
!密度：噸/立方毫米 能量：牛頓*毫米=0。001焦耳； 加速度：MM/S^-2=0.001g.

!****************************************************************************

finish
/clear,nostart
/filename,yhbridge,1
/title,yhbridge
/prep7

!NSKIP里存儲的是本次調(diào)用在迭代序列中的序號，ANSYS朝FREQUENCE里保存數(shù)據(jù)時才用追加的形式保存
!而MATLAB朝OPTPARA里保存數(shù)據(jù)也是采用追加的方式保存。
!這樣，當優(yōu)化運行結(jié)束時，OPTPARA 里每行儲存的數(shù)據(jù)就對應設計空間內(nèi)一個相應的點，
!而FREQUENCE 則是該點對應的各階頻率 。

! NSLOOP 里按照逐行順序，存儲的是該次調(diào)用時各個循環(huán)指標的值（每個循環(huán)指標對應一個修正系數(shù)OPT）。
! NSLOOP 是按照刷新的方式存儲的，所以里面的數(shù)據(jù)僅僅對應一個設計點（對應最后一次調(diào)用）。

!而NSKIP則是該次循環(huán)的序號，顯然是各個指標的乘積。
*dim,NSLOOP,table,3,4
*tread,NSLOOP,NSLOOP,dat,'e:\ansyswork\yhdata\input'
NSKIP=1
*do,i,1,3,1
*do,j,1,4,1
NSKIP=NSKIP*NSLOOP(i,j)
*enddo
*enddo

! 由于是采用非追加的方式保存數(shù)據(jù)，所以，直接讀出數(shù)據(jù)即可。
*dim,optpara,table,1,12,1
*tread,optpara,optpara,dat,'e:\ansyswork\yhdata\input',0
! opt1修正主梁彈性模量，opt2 修正主梁密度；opt3修正橋塔彈性模量，opt4修正橋塔密度；opt5修正索的彈性模量.
! pt6修正主塔垂直橋縱軸向水平約束；opt7 對應橋頭水平約束的修正系數(shù)；opt8 對應輔助橋墩的水平約束修正系數(shù)。
!pt9 修正橋縱軸向水平約束；pt10 修正抗扭約束；pt11 修正4個橋墩豎向約束，pt12修正索幅在橋面上的寬度（橋面位置上下游索間距）

opt1=optpara(1,1)
opt2=optpara(1,2)
opt3=optpara(1,3)
opt4=optpara(1,4)
opt5=optpara(1,5)
opt6=optpara(1,6)
opt7=optpara(1,7)
opt8=optpara(1,8)
opt9=optpara(1,9)
opt10=optpara(1,10)
opt11=optpara(1,11)
opt12=optpara(1,12)

!定義主梁為第一類單元、實常數(shù)、材料組
et,1,189
mp,prxy,1,0.17
mp,ex,1,opt1*3.5e4
mp,dens,1,opt2*1.01*2.4e-9

!定義魚刺（橫向梁）為第二類單元 、實常數(shù)、材料組
et,2,4
mp,prxy,2,0.01
mp,ex,2,3.5e4
dens,2,0.002e-9

!定義橋塔為第三類單元、實常數(shù)、材料組
et,3,188
mp,prxy,3,0.17
mp,ex,3,opt3*0.95*3.5e4
dens,3,opt4*2.4e-9

!定義索為第四類單元。由于實常數(shù)、及彈性模量各個索有區(qū)別，所以在后面分別定義。
et,4,8

!定義主梁的彈性水平約束，采用 LINK 8 單元等效彈性約束
et,5,8
mp,prxy,5,0.3
mp,ex,5,2.0e5
dens,5,0.002e-9

!定義所有beam 4 單元的實常數(shù)。 注意實常數(shù)與單位有關(guān)。剛臂的實常數(shù)為R2。主梁當采用BEAM 4
!單元時，實在常數(shù)為R1
r,1,0.684e7,0.147e15,0.221e13,100,12000,,,,,,
r,2,1.0e6,0.33e14,0.33e14,100,100

!r,3,0.684e7,0,
!3 號實常數(shù)為主橋塔的垂直橋縱軸的水平彈性約束等效成拉壓桿件的面積 。
!4 號實常數(shù)為兩側(cè)橋頭的垂直橋縱軸的水平彈性約束等效成拉壓桿件的面積 。
!5 號實常數(shù)為兩側(cè)輔助橋墩的垂直橋縱軸的水平彈性約束等效成拉壓桿件的面積 。
r,3,opt6*52.5,0
r,4,opt7*1100,0
r,5,opt8*1100,0

!6 號實常數(shù)為兩側(cè)橋頭的沿橋縱軸向彈性約束等效成拉壓桿件的面積 。

r,6,opt9*1,0
r7=opt10*40000
r8=opt11*90000-2*r7

!7 號實常數(shù)為兩側(cè)橋頭的外側(cè)豎向抗扭轉(zhuǎn)彈性約束等效成拉壓桿件的面積 。
!7 號實常數(shù)為兩側(cè)輔助橋墩的外側(cè)抗扭轉(zhuǎn)彈性約束等效成拉壓桿件的面積 。
r,7,r7,0
!8 號實常數(shù)為兩側(cè)輔助橋墩及橋頭墩的中心豎向彈性約束等效成拉壓桿件的面積 。
r,8,r8,0

!*************** 下面定義索的實常數(shù) ***********************

!從CABLEDATA 里讀出索的面積及初始應變。
!其中第一列為河側(cè)索面積，單位為平方毫米。第二列為河側(cè)索初始索力，單位為噸。第三列是岸側(cè)索面積
!第四列是岸側(cè)索初始索力。第五列是河側(cè)索的水平投影，第六列是岸側(cè)索的水平投影.當考慮索的非線性，
!采用 ERNST 修正法來修正索的彈性模量，要注意單位轉(zhuǎn)化。設置索的容重gama，來作為后面計算索的彈性
!模量的數(shù)據(jù)。容重單位為KN/M^3,彈性模量單位為帕。
!定義數(shù)組 Eeqs 來存儲后面生成的等效彈性模量。 第一列存儲河側(cè)索的等效彈性模量,第二列存儲岸側(cè)索的
!等效彈性模量。定義 COFF 來存儲彈性模量等效時的系數(shù)。這樣更容易看出是否出錯。這個系數(shù)應該是略
!微小于1的正數(shù)。
!現(xiàn)在索的初始張拉力已經(jīng)全部在后面用 CABLEF 來存儲了，CABLEDATA 里的初始索力數(shù)據(jù)為**研究
!所的設計數(shù)據(jù).

*dim,rcable,table,11,6,1
*tread,rcable,cabledata,txt,'e:\ansyswork\yhdata\input',

*dim,cablef,table,11,8,1
*tread,cablef,cablef,txt,'e:\ansyswork\yhdata\input',1

*set,gama,78000
*SET,e0,2.07e11

*dim,Eeqs,array,11,8,1
*dim,coff,array,11,8,1

!為了后面模擬施工過程的換索方便，將每根索都單獨賦實常數(shù)。
!為了模擬有無初始索力的方便，在索的初始應變前都有個系數(shù)BEITA，
!當不考慮初始應變的時候，就給BEITA 賦 0 值。考慮，則為 1 。
!上游部分：
!定義北塔岸側(cè)索的材料及實常數(shù) 。 編號從 11 開始到 21 結(jié)束
!定義北塔河側(cè)索的材料及實常數(shù) 。 編號從 31 開始到 41 結(jié)束
!定義北塔河側(cè)索的材料及實常數(shù) 。 編號從 51 開始到 61 結(jié)束
!定義南塔岸側(cè)索的材料及實常數(shù) 。 編號從 71 開始到 81 結(jié)束
!下游部分：
!定義北塔岸側(cè)索的材料及實常數(shù) 。 編號從 91 開始到 101 結(jié)束
!定義北塔河側(cè)索的材料及實常數(shù) 。 編號從 111 開始到 121 結(jié)束
!定義北塔河側(cè)索的材料及實常數(shù) 。 編號從 131 開始到 141 結(jié)束
!定義南塔岸側(cè)索的材料及實常數(shù) 。 編號從 151 開始到 161 結(jié)束

*do,i,1,11,1
*set,coff1,(rcable(i,6)*0.001*gama)**2
*set,coff2,12*(cablef(i,1)*1.0e6/rcable(i,3))**3
coff(i,1)=1/((coff1*e0/coff2)+1)
Eeqs(i,1)=opt5*0.79*coff(i,1)*e0*1.0e-6
mp,prxy,i+10,0.3
mp,ex,i+10,Eeqs(i,1)
mp,dens,i+10,0.1e-11
mp,alpx,i+10,0.00001
*enddo

*do,i,1,11,1
*set,coff1,(rcable(i,5)*0.001*gama)**2
*set,coff2,12*(cablef(i,2)*1.0e6/rcable(i,1))**3
coff(i,2)=1/((coff1*e0/coff2)+1)
Eeqs(i,2)=opt5*0.79*coff(i,2)*e0*1.0e-6
mp,prxy,i+30,0.3
mp,ex,i+30,Eeqs(i,2)
mp,dens,i+30,0.1e-11
mp,alpx,i+30,0.00001
*enddo

*do,i,1,11,1
*set,coff1,(rcable(i,5)*0.001*gama)**2
*set,coff2,12*(cablef(i,3)*1.0e6/rcable(i,1))**3
coff(i,3)=1/((coff1*e0/coff2)+1)
Eeqs(i,3)=opt5*0.79*coff(i,3)*e0*1.0e-6
mp,prxy,i+50,0.3
mp,ex,i+50,Eeqs(i,3)
mp,dens,i+50,0.1e-11
mp,alpx,i+50,0.00001
*enddo

*do,i,1,11,1
*set,coff1,(rcable(i,6)*0.001*gama)**2
*set,coff2,12*(cablef(i,4)*1.0e6/rcable(i,3))**3
coff(i,4)=1/((coff1*e0/coff2)+1)
Eeqs(i,4)=opt5*0.79*coff(i,4)*e0*1.0e-6
mp,prxy,i+70,0.3
mp,ex,i+70,Eeqs(i,4)
mp,dens,i+70,0.1e-11
mp,alpx,i+70,0.00001
*enddo

*do,i,1,11,1
*set,coff1,(rcable(i,6)*0.001*gama)**2
*set,coff2,12*(cablef(i,5)*1.0e6/rcable(i,3))**3
coff(i,5)=1/((coff1*e0/coff2)+1)
Eeqs(i,5)=opt5*0.79*coff(i,5)*e0*1.0e-6
mp,prxy,i+90,0.3
mp,ex,i+90,Eeqs(i,5)
mp,dens,i+90,0.1e-11
mp,alpx,i+90,0.00001
*enddo

*do,i,1,11,1
*set,coff1,(rcable(i,5)*0.001*gama)**2
*set,coff2,12*(cablef(i,6)*1.0e6/rcable(i,1))**3
coff(i,6)=1/((coff1*e0/coff2)+1)
Eeqs(i,6)=opt5*0.79*coff(i,6)*e0*1.0e-6
mp,prxy,i+110,0.3
mp,ex,i+110,Eeqs(i,6)
mp,dens,i+110,0.1e-11
mp,alpx,i+110,0.00001
*enddo

*do,i,1,11,1
*set,coff1,(rcable(i,5)*0.001*gama)**2
*set,coff2,12*(cablef(i,7)*1.0e6/rcable(i,1))**3
coff(i,7)=1/((coff1*e0/coff2)+1)
Eeqs(i,7)=opt5*0.79*coff(i,7)*e0*1.0e-6
mp,prxy,i+130,0.3
mp,ex,i+130,Eeqs(i,7)
mp,dens,i+130,0.1e-11
mp,alpx,i+130,0.00001
*enddo

*do,i,1,11,1
*set,coff1,(rcable(i,6)*0.001*gama)**2
*set,coff2,12*(cablef(i,8)*1.0e6/rcable(i,3))**3
coff(i,8)=1/((coff1*e0/coff2)+1)
Eeqs(i,8)=opt5*0.79*coff(i,8)*e0*1.0e-6
mp,prxy,i+150,0.3
mp,ex,i+150,Eeqs(i,8)
mp,dens,i+150,0.1e-11
mp,alpx,i+150,0.00001
*enddo

!為了后面模擬施工過程的換索方便，將每根索都單獨賦實常數(shù)。
!為了模擬有無初始索力的方便，在索的初始應變前都有個系數(shù)BEITA，
!當不考慮初始應變的時候，就給BEITA 賦 0 值。考慮，則為 1 。
!再賦實常數(shù)。
beita=1
*do,i,1,11,1
R,i+10,rcable(i,3),beita*cablef(i,1)/rcable(i,3)/Eeqs(i,1)
R,i+30,rcable(i,1),beita*cablef(i,2)/rcable(i,1)/Eeqs(i,2)
R,i+50,rcable(i,1),beita*cablef(i,3)/rcable(i,1)/Eeqs(i,3)
R,i+70,rcable(i,3),beita*cablef(i,4)/rcable(i,3)/Eeqs(i,4)
R,i+90,rcable(i,3),beita*cablef(i,5)/rcable(i,3)/Eeqs(i,5)
R,i+110,rcable(i,1),beita*cablef(i,6)/rcable(i,1)/Eeqs(i,6)
R,i+130,rcable(i,1),beita*cablef(i,7)/rcable(i,1)/Eeqs(i,7)
R,i+150,rcable(i,3),beita*cablef(i,8)/rcable(i,3)/Eeqs(i,8)
*enddo

!下面開始定義塔的幾何模型。坐標起點從對稱縱剖面、橫剖面與基礎(chǔ)承臺的交點為坐標原點。
!下面定義塔腳坐標
k,1,8620,5000
k,2,-8620,5000

!下面定義塔第一橫梁中心坐標。
k,3,6250,44758
k,4,-6250,44758
!下面定義兩側(cè)雙排索從第一個索到第11個索的坐標。
k,5,6250,46100
k,6,-6250,46100

k,7,6250,48850
k,8,-6250,48850
k,9,6250,50450
k,10,-6250,50450
k,11,6250,51750
k,12,-6250,51750
k,13,6250,53000
k,14,-6250,53000
k,15,6250,54200
k,16,-6250,54200
k,17,6250,55400
k,18,-6250,55400
k,19,6250,56600
k,20,-6250,56600
k,21,6250,57800
k,22,-6250,57800
k,23,6250,59000
k,24,-6250,59000
!下面定義塔第二個橫梁的坐標。
k,25,6250,59625
k,26,-6250,59625
!下面定義塔頂坐標。
k,27,6250,65000
k,28,-6250,65000

!下面定義塔截面的方向關(guān)鍵點。
k,31,6250,44758,-10000
k,32,-6250,44758,-10000
k,33,6250,44758,-10000
k,34,6250,59625,-10000

!由于1、2兩根短索接近垂直，與塔的夾角很小。直接將其朝軸線上簡化誤差比較大。
!為次，在該截面建立兩個剛臂來解決。下面定義剛臂端點節(jié)點的關(guān)鍵點。

!北塔上游的剛臂
k,50,-6250,46000,126200-1200
k,51,-6250,46000,126200+1200
!北塔下游的剛臂
k,52,6250,46000,126200-1200
k,53,6250,46000,126200+1200
!南塔上游的剛臂
k,54,-6250,46000,386200-1200
k,55,-6250,46000,386200+1200
!南塔下游的剛臂
k,56,6250,46000,386200-1200
k,57,6250,46000,386200+1200

cm,ktower,kp
cmsel,s,ktower
!復制生成北橋塔的關(guān)鍵點，編號增加量為100，以便于控制。
kgen,2,1,34,1,0,0,126200,100,0,0
!移動生成南橋塔，編號增加量再增加100，以便于控制。
kgen,2,1,34,1,0,0,386200,0,0,1

!下面將塔的關(guān)鍵點用線連接起來，生成橋梁主塔的中心線。

*do,i,1,25,2
l,i,i+2
l,i+1,i+3
*enddo
l,3,4
l,25,26

numstr,line,101
*do,i,101,125,2
l,i,i+2
l,i+1,i+3
*enddo
l,103,104
l,125,126

!整個南塔北塔已經(jīng)建模完畢。現(xiàn)在開始建立主梁模型。主糧模型按照單梁魚脊骨模型
!建立。整個梁自然節(jié)點的Z坐標已經(jīng)保存在數(shù)據(jù)文件 beamdata.txt 里面了。該數(shù)據(jù)的格式為：
!第一行和第一列分別表示行標號、列標號；按照列順序分別為X、Y、Z坐標。讀入該數(shù)據(jù)到
!數(shù)組參數(shù) axis ，再用beamb、bdeep 兩個標量參數(shù)來表示等效的索橫向間距（橋面計算寬度）
!和斜拉索偏離主梁形心的間距，以考慮斜拉索建立的結(jié)構(gòu)預應力對主梁產(chǎn)生的反拱。

*set,beamb,12000*opt12
*set,bdeep,368.25

!定義 axis 數(shù)組的維數(shù)。由于主梁在該計算模型中被分為51個節(jié)點，所以該數(shù)組
!為51*3的數(shù)組。

*dim,axis,table,51,3,1
*tread,axis,beamdata,txt,'e:\ansyswork\yhdata\input',

!以上 axis 里的是橋面的坐標及高程數(shù)據(jù)，但是實際上主梁的中性軸線距離橋面距離為613.754毫米，所以
!還需要把 axis 的第二列，也就是Y坐標向下平移這個距離，使得 axis 里是縱軸所在的坐標。

*do,i,1,11
axis(i,2)=axis(i,2)-613.754
*enddo

!定義cabup 數(shù)組的維數(shù)。該數(shù)組保存的是上游索的下節(jié)點位置。只保存Z坐標及Y坐標。X坐標由半幅梁寬得到。
!為51*3的數(shù)組。
*dim,cabup,table,51,3,1
cabup(0,0)=0
cabup(0,1)=1
cabup(0,2)=2
cabup(0,3)=3
*do,i,1,51,1
cabup(i,0)=i
cabup(i,1)=-1*beamb/2
cabup(i,2)=axis(i,2)-bdeep
cabup(i,3)=axis(i,3)
*enddo

!定義cabdown 數(shù)組的維數(shù)。該數(shù)組保存的是下游索的下節(jié)點位置。 只保存Z坐標及Y坐標。X坐標由半幅梁寬得到。
*dim,cabdown,table,51,3,1
cabdown(0,0)=0
cabdown(0,1)=1
cabdown(0,2)=2
cabdown(0,3)=3
*do,i,1,51,1
cabdown(i,0)=i
cabdown(i,1)=1*beamb/2
cabdown(i,2)=axis(i,2)-bdeep
cabdown(i,3)=axis(i,3)
*enddo
!定義bnup 數(shù)組的維數(shù)。該數(shù)組保存的是魚脊椎上游方向節(jié)點位置。 只保存Z坐標及Y坐標。X坐標由半幅梁寬得到。
*dim,bnup,table,51,3,1
bnup(0,0)=0
bnup(0,1)=1
bnup(0,2)=2
bnup(0,3)=3
*do,i,1,51,1
bnup(i,0)=i
bnup(i,1)=-1*beamb/2
bnup(i,2)=axis(i,2)
bnup(i,3)=axis(i,3)
*enddo
!定義bndown 數(shù)組的維數(shù)。該數(shù)組保存的是魚脊椎上游方向節(jié)點位置。 只保存Z坐標及Y坐標。X坐標由半幅梁寬得到。
*dim,bndown,table,51,3,1
bndown(0,0)=0
bndown(0,1)=1
bndown(0,2)=2
bndown(0,3)=3
*do,i,1,51,1
bndown(i,0)=i
bndown(i,1)=1*beamb/2
bndown(i,2)=axis(i,2)
bndown(i,3)=axis(i,3)
*enddo
!上面已經(jīng)建立了梁的全部關(guān)鍵節(jié)點坐標。下面將利用上面各個參數(shù)數(shù)組里保存的
!相關(guān)坐標數(shù)據(jù)來生成魚脊椎梁模型的線。生成的節(jié)點將從301開始編號。 300號節(jié)點
!將用來作為主梁的截面方向關(guān)鍵點。

!生成梁的軸線上關(guān)鍵節(jié)點。方向關(guān)鍵點K300。

k,300,0,30000,20000
*do,i,1,51,1
k,i+300,axis(i,1),axis(i,2),axis(i,3)
*enddo
!***********************************************************************
!此處采用 LINK 8 單元 生成所有構(gòu)成邊界條件的彈性約束。

!生成梁在橋墩處的垂直軸線方向的水平彈性約束，采用 LINK8 單元。
!spring 是邊界的彈性桿長度，為了美觀可以調(diào)節(jié)該參數(shù)。
spring=1000
! 355，356 ，357，358四個點分別在四個橋墩位置的正下方，用來生成豎向約束（模擬盆式支座）
k,355,axis(1,1),axis(1,2)-spring,axis(1,3)
k,356,axis(5,1),axis(5,2)-spring,axis(5,3)
k,357,axis(47,1),axis(47,2)-spring,axis(47,3)
k,358,axis(51,1),axis(51,2)-spring,axis(51,3)
! 361，362 ，363，364四個點分別在四個橋墩位置的側(cè)向，用來生成水平約束（模擬盆式支座的摩擦力約束）
k,361,axis(1,1)+spring,axis(1,2),axis(1,3)
k,362,axis(5,1)+spring,axis(5,2),axis(5,3)
k,363,axis(47,1)+spring,axis(47,2),axis(47,3)
k,364,axis(51,1)+spring,axis(51,2),axis(51,3)

!生成梁在橋墩處的沿橋軸線方向的水平彈性約束，采用 LINK8 單元，用來生成水平約束（模擬盆式支座的摩擦力約束）。
k,371,axis(1,1),axis(1,2),axis(1,3)-spring
!生成梁在主塔處的沿橋垂直軸線方向的水平彈性約束，采用 LINK8 單元，用來生成水平約束（模擬橋塔的碰撞）。
k,372,axis(14,1)+spring,axis(14,2),axis(14,3)
k,373,axis(38,1)+spring,axis(38,2),axis(38,3)

!梁在橋墩處邊上的的垂直方向的抗扭轉(zhuǎn)彈性約束，采用 LINK8 單元。但是不需要再次建立線
!因為后面建立索端剛臂時在編號含有1，5，47，51 處時，恰好不是索端而是橋墩位置，在后
!面劃分網(wǎng)格的時候注意把這些打算作為邊界條件彈性約束的線挑選出來劃分成邊界條件的單元
!而非索端剛臂！
!***********************************************************************

! 生成梁上游的全部脊椎骨端節(jié)點 。節(jié)點編號從401開始
*do,i,1,51,1
k,i+400,bnup(i,1),bnup(i,2),bnup(i,3)
*enddo
!生成梁下游的全部脊椎骨端節(jié)點 。節(jié)點編號從501開始
*do,i,1,51,1
k,i+500,bndown(i,1),bndown(i,2),bndown(i,3)
*enddo

! 生成梁上游的全部索端節(jié)點 。節(jié)點編號從601開始 ，

k,1+600,cabup(1,1),cabup(1,2)-spring,cabup(1,3)

*do,i,2,4,1
k,i+600,cabup(i,1),cabup(i,2),cabup(i,3)
*enddo

k,5+600,cabup(5,1),cabup(5,2)-spring,cabup(5,3)

*do,i,6,46,1
k,i+600,cabup(i,1),cabup(i,2),cabup(i,3)
*enddo

k,47+600,cabup(47,1),cabup(47,2)-spring,cabup(47,3)

*do,i,48,50,1
k,i+600,cabup(i,1),cabup(i,2),cabup(i,3)
*enddo

k,51+600,cabup(51,1),cabup(51,2)-spring,cabup(51,3)

!生成梁下游的全部索端節(jié)點 。節(jié)點編號從701開始 ，

k,1+700,cabdown(1,1),cabdown(1,2)-spring,cabdown(1,3)

*do,i,2,4,1
k,i+700,cabdown(i,1),cabdown(i,2),cabdown(i,3)
*enddo

k,5+700,cabdown(5,1),cabdown(5,2)-spring,cabdown(5,3)

*do,i,6,46,1
k,i+700,cabdown(i,1),cabdown(i,2),cabdown(i,3)
*enddo

k,47+700,cabdown(47,1),cabdown(47,2)-spring,cabdown(47,3)

*do,i,48,50,1
k,i+700,cabdown(i,1),cabdown(i,2),cabdown(i,3)
*enddo

k,51+700,cabdown(51,1),cabdown(51,2)-spring,cabdown(51,3)

!生成梁的軸線 。編號從300開始編號。 編號為從300到349號。

numstr,line,300
*do,i,1,50,1
l,i+300,i+301
*enddo

!由前面描述的設計數(shù)據(jù)，301，305，347，351 分別是橋主梁軸線在橋墩位置處的坐標。
numstr,line,361
l,361,301
l,362,305
l,363,347
l,364,351

!生成橋縱軸向水平彈性約束，編號為370號。
numstr,line,370
l,301,371
l,314,372
l,338,373

!生成橋豎向水平彈性約束，編號為380號。
numstr,line,380
l,301,355
l,305,356
l,347,357
l,351,358

!生成全部脊椎骨上的魚刺梁。 編號從400開始編號。
numstr,line,400
*do,i,1,51,1
l,i+300,i+400
l,i+300,i+500
*enddo

!生成全部上游索端剛臂梁。 編號從600開始編號。
numstr,line,600
*do,i,1,51,1
l,i+400,i+600
*enddo

!生成全部下游索端剛臂梁。 編號從700開始編號。
numstr,line,700
*do,i,1,51,1
l,i+500,i+700
*enddo

!以上已經(jīng)全部建立完畢梁單元和塔單元。下面將建立索單元。

!主梁北塔北側(cè)上游索端點編號我們從最短索到最長索存到 BCABLE1 。
!主梁北塔南側(cè)上游索端點編號我們從最短索到最長索存到 BCABLE2 。
!主梁南塔北側(cè)上游索端點編號我們從最短索到最長索存到 BCABLE3 。
!主梁南塔南側(cè)上游索端點編號我們從最短索到最長索存到 BCABLE4 。

!主梁北塔北側(cè)下游索端點編號我們從最短索到最長索存到 BCABLE5 。
!主梁北塔南側(cè)下游索端點編號我們從最短索到最長索存到 BCABLE6 。
!主梁南塔北側(cè)下游索端點編號我們從最短索到最長索存到 BCABLE7 。
!主梁南塔南側(cè)下游索端點編號我們從最短索到最長索存到 BCABLE8 。

!下面開始定義以上12個參數(shù)的維數(shù)。
*dim,BCABLE1,array,11,1,1
*dim,BCABLE2,array,11,1,1
*dim,BCABLE3,array,11,1,1
*dim,BCABLE4,array,11,1,1
*dim,BCABLE5,array,11,1,1
*dim,BCABLE6,array,11,1,1
*dim,BCABLE7,array,11,1,1
*dim,BCABLE8,array,11,1,1

!邊界條件的說明：
!北塔北側(cè) 1 號橋墩的節(jié)點號碼為1（301、401、501、601、701）
!北塔北側(cè) 2 號橋墩的節(jié)點號碼為5（305、405、505、605、705）
!北塔 3 號橋墩（主塔自身）的節(jié)點號碼為14（314、414、514、614、714）

!中心對稱點的節(jié)點號碼為 26 。
!北塔 4 號橋墩（主塔自身）的節(jié)點號碼為38 。
!南塔南側(cè) 5 號橋墩的節(jié)點號碼為47 。
!南塔南側(cè) 5 號橋墩的節(jié)點號碼為51 。

BCABLE1(1)=613,612,611,610,609,608,607,606
BCABLE1(9)=604,603,602

BCABLE2(1)=615,616,617,618,619,620,621,622
BCABLE2(9)=623,624,625

BCABLE3(1)=637,636,635,634,633,632,631,630
BCABLE3(9)=629,628,627

BCABLE4(1)=639,640,641,642,643,644,645,646
BCABLE4(9)=648,649,650
*do,i,1,11,1
BCABLE5(i,1)=BCABLE1(i,1)+100
BCABLE6(i,1)=BCABLE2(i,1)+100
BCABLE7(i,1)=BCABLE3(i,1)+100
BCABLE8(i,1)=BCABLE4(i,1)+100
*enddo

!下面再來使用查詢函數(shù)KP來得到橋塔上索節(jié)點的編號，并把他們保存到TCABLE 數(shù)組總?cè)ァ?br>!北塔上游塔上索節(jié)點從最短索到最長索的編號我們將收集來存在數(shù)組 TCABLE1 。
!北塔下游塔上索節(jié)點從最短索到最長索的編號我們將收集來存在數(shù)組 TCABLE2 。
!南塔上游塔上索節(jié)點從最短索到最長索的編號我們將收集來存在數(shù)組 TCABLE3 。
!南塔下游塔上索節(jié)點從最短索到最長索的編號我們將收集來存在數(shù)組 TCABLE4 。

*dim,tcable1,array,11,1,1
*dim,tcable2,array,11,1,1
*dim,tcable3,array,11,1,1
*dim,tcable4,array,11,1,1
*dim,cabcoord,array,11,1,1

! 北橋塔下游索中心點的坐標為（6250，，126200）
! 北橋塔上游索中心點的坐標為（-6250，，126200）
! 南橋塔索下游中心點的坐標為（6250，，386200）
! 南橋塔索上游中心點的坐標為（-6250，，386200）
! 南橋塔索中心點的坐標見下列數(shù)組 cabcoord 中。

cabcoord(1)=59000,57800,56600,55400,54200,53000,51750,50450
cabcoord(9)=48850,46000,46000

*do,i,1,11,1
num=12-i
tcable1(num)=kp(-62500,cabcoord(i),126200)
tcable2(num)=kp(62500,cabcoord(i),126200)
tcable3(num)=kp(-62500,cabcoord(i),386200)
tcable4(num)=kp(62500,cabcoord(i),386200)
*enddo

!生成全部塔上連接1、2兩號索的小剛臂梁。 編號從801開始。
numstr,line,801
l,50,tcable1(1)
l,51, tcable1(1)
l,52,tcable2(1)
l,53,tcable2(1)
l,54,tcable3(1)
l,55,tcable3(1)
l,56,tcable4(1)
l,57,tcable4(1)

!下面用上面 TCABLE 、BCABLE 里面的數(shù)據(jù)來建立索單元。索單元的編號規(guī)則如下：
!從1001號開始，每排索面為一組，
!第一組：北塔北上游 。號碼為1001-1011
!第二組：北塔南上游 。號碼為1101-1111
!第三組：南塔北上游 。號碼為1201-1211
!第四組：南塔南上游 。號碼為1301-1311
!第五組：北塔北下游 。號碼為1401-1411
!第六組：北塔南下游 。號碼為1501-1511
!第七組：南塔北下游 。號碼為1601-1611
!第八組：南塔南下游 。號碼為1701-1711

!建立第一組索單元。北塔北上游
numstr,line,1001
l,50,bcable1(1)
l,50,bcable1(2)
*do,i,3,11,1
l,tcable1(i),bcable1(i)
*enddo
!建立第二組索單元 。北塔南上游 。
numstr,line,1101
l,51,bcable2(1)
l,51,bcable2(2)
*do,i,3,11,1
l,tcable1(i),bcable2(i)
*enddo
!建立第三組索單元 。南塔北上游。
numstr,line,1201
l,54,bcable3(1)
l,54,bcable3(2)
*do,i,3,11,1
l,tcable3(i),bcable3(i)
*enddo
!建立第四組索單元。 南塔南上游。
numstr,line,1301
l,55,bcable4(1)
l,55,bcable4(2)
*do,i,3,11,1
l,tcable3(i),bcable4(i)
*enddo
!建立第五組索單元 。北塔北下游 。
numstr,line,1401
l,52,bcable5(1)
l,52,bcable5(2)
*do,i,3,11,1
l,tcable2(i),bcable5(i)
*enddo
!建立第六組索單元 。北塔南下游 。
numstr,line,1501
l,53,bcable6(1)
l,53,bcable6(2)
*do,i,3,11,1
l,tcable2(i),bcable6(i)
*enddo
!建立第七組索單元 。南塔北下游 。
numstr,line,1601
l,56,bcable7(1)
l,56,bcable7(2)
*do,i,3,11,1
l,tcable4(i),bcable7(i)
*enddo
!建立第八組索單元 。南塔南下游 。
numstr,line,1701
l,57,bcable8(1)
l,57,bcable8(2)
*do,i,3,11,1
l,tcable4(i),bcable8(i)
*enddo
SAVE
!整個前處理單元建立完畢

!************************************************************************************************

!永和橋網(wǎng)格劃分

!************************************************************************************************
!下面將塔的截面全部讀進內(nèi)存。
sectype,1,beam,mesh !指定讀入的截面類型在后面使用中編號
secoffset,shrc !指定截面在梁縱軸上的偏移量
secread,'yhtsec1','sect','e:\ansyswork\yhbridge\seclib',mesh !讀入截面。如果截面保存在其他路徑，可以采用絕對路徑的方法確定

sectype,2,beam,mesh !指定讀入的截面類型在后面使用中編號
secoffset,shrc!指定截面在梁縱軸上的偏移量
secread,'yhtsec2','sect','e:\ansyswork\yhbridge\seclib',mesh !讀入截面。如果截面保存在其他路徑，可以采用絕對路徑的方法確定

sectype,3,beam,mesh !指定讀入的截面類型在后面使用中編號
secoffset,shrc !指定截面在梁縱軸上的偏移量
secread,'yhtsec3','sect','e:\ansyswork\yhbridge\seclib',mesh !讀入截面。如果截面保存在其他路徑，可以采用絕對路徑的方法確定

sectype,4,beam,mesh !指定讀入的截面類型在后面使用中編號
secoffset,shrc !指定截面在梁縱軸上的偏移量
secread,'yhtsec4','sect','e:\ansyswork\yhbridge\seclib',mesh !讀入截面。如果截面保存在其他路徑，可以采用絕對路徑的方法確定

sectype,5,beam,mesh !指定讀入的截面類型在后面使用中編號
secoffset,cent !指定截面在梁縱軸上的偏移量
secread,'yhbsec1','sect','e:\ansyswork\yhbridge\seclib',mesh !讀入截面。如果截面保存在其他路徑，可以采用絕對路徑的方法確定

!對于188\189單元，不需要實常數(shù)，但是如果不給其指定一個空的實常數(shù)號，會自動把第一個實常數(shù)號給它分配。產(chǎn)生警告信息。

!***********索的分布編號說明*************
!北塔岸側(cè)索：1001-1011；1401-1411
!北塔河側(cè)索：1101-1111；1501-1511
!南塔河側(cè)索：1201-1211；1601-1411
!南塔岸側(cè)索：1301-1311；1701-1711
!按照從上游到下游，從北到南的順序！

numstr,elem,1
*do,i,1,11
lsel,s,line,,1000+i
latt,i+10,i+10,4
esize,,1
lmesh,all
*enddo

numstr,elem,101
*do,i,1,11
lsel,s,line,,1100+i
latt,i+30,i+30,4
esize,,1
lmesh,all
*enddo

numstr,elem,201
*do,i,1,11
lsel,s,line,,1200+i
latt,i+50,i+50,4
esize,,1
lmesh,all
*enddo

numstr,elem,301
*do,i,1,11
lsel,s,line,,1300+i
latt,i+70,i+70,4
esize,,1
lmesh,all
*enddo

numstr,elem,401
*do,i,1,11
lsel,s,line,,1400+i
latt,i+90,i+90,4
esize,,1
lmesh,all
*enddo

numstr,elem,501
*do,i,1,11
lsel,s,line,,1500+i
latt,i+110,i+110,4
esize,,1
lmesh,all
*enddo

numstr,elem,601
*do,i,1,11
lsel,s,line,,1600+i
latt,i+130,i+130,4
esize,,1
lmesh,all
*enddo

numstr,elem,701
*do,i,1,11
lsel,s,line,,1700+i
latt,i+150,i+150,4
esize,,1
lmesh,all
*enddo

! 再給四個橋墩位置的橫橋向水平彈性約束劃分單元，同樣采用LINK 8 單元。
! 采用編號控制，其LINE 的編號為361-364

numstr,elem,801

lsel,s,line,,361,361,1,0
latt,5,4,5
esize,,1
lmesh,all

lsel,s,line,,362,362,1,0
latt,5,5,5
esize,,1
lmesh,all
lsel,s,line,,363,363,1,0
latt,5,5,5
esize,,1
lmesh,all
lsel,s,line,,364,364,1,0
latt,5,4,5
esize,,1
lmesh,all

! 再給橋墩位置的橋縱軸向水平彈性約束劃分單元，同樣采用LINK 8 單元。
! 采用編號控制，其LINE 的編號為370

lsel,s,line,,370,370,1,0
latt,5,6,5
esize,,1
lmesh,all

! 再給住塔位置的橋垂直軸向水平彈性約束劃分單元，同樣采用LINK 8 單元。
! 采用編號控制，其LINE 的編號為371,372

lsel,s,line,,371,372,1,0
latt,5,3,5
esize,,1
lmesh,all

numstr,elem,1000
!開始給魚脊椎骨主梁劃分單元
lsel,s,line,,300,349,1,0
latt,1,1,1,,300,,5
esize,2000
lmesh,all

numstr,elem,2000
!開始給魚脊剛性橫梁主梁劃分單元 .
lsel,s,line,,400,501,1,0
latt,2,2,2,,,,
esize,,1
lmesh,all

!開始給魚脊椎骨下的剛性剛臂梁劃分單元 .
lsel,s,line,,601,649,1,0
lsel,a,line,,701,749,1,0
lsel,u,line,,704,704,1,0
lsel,u,line,,604,604,1,0
lsel,u,line,,646,646,1,0
lsel,u,line,,746,746,1,0
latt,2,2,2,,,,
esize,,1
lmesh,all

!開始給邊界條件的扭轉(zhuǎn)彈簧劃分單元中間兩個輔助墩.
lsel,s,line,,704,704,1,0
lsel,a,line,,604,604,1,0
lsel,a,line,,646,646,1,0
lsel,a,line,,746,746,1,0
latt,5,7,5,,,,
esize,,1
lmesh,all

!開始給邊界條件的扭轉(zhuǎn)彈簧劃分單元。 兩個橋頭橋墩.
lsel,s,line,,700,700,1,0
lsel,a,line,,600,600,1,0
lsel,a,line,,650,650,1,0
lsel,a,line,,750,750,1,0
latt,5,7,5,,,,
esize,,1
lmesh,all

!開始給邊界條件的豎向彈簧劃分單元。
lsel,s,line,,380,380,1,0
lsel,a,line,,381,381,1,0
lsel,a,line,,382,382,1,0
lsel,a,line,,383,383,1,0
latt,5,8,5,,,,
esize,,1
lmesh,all

!開始給塔上的小剛臂梁劃分單元 。
lsel,s,line,,801,808,1,0
latt,2,2,2,,,,
esize,,1
lmesh,all

numstr,elem,4000

!先給北塔分配屬性和單元劃分，區(qū)別僅僅是關(guān)鍵點及線的編號增量為100。
lsel,s,line,,101,101,1,0
latt,3,25,3,,131,,1
esize,2000
lmesh,all

lsel,s,line,,102,102,1,0
latt,3,25,3,,132,,1
esize,2000
lmesh,all

lsel,s,line,,103,125,2,0
latt,3,25,3,,131,,3
esize,2000
lmesh,all

lsel,s,line,,104,126,2,0
latt,3,25,3,,132,,3
esize,2000
lmesh,all

lsel,s,line,,127,127,1,0
latt,3,25,3,,133,,4
esize,2000
lmesh,all

lsel,s,line,,128,128,1,0
latt,3,25,3,,134,,4
esize,2000
lmesh,all

numstr,elem,5000
!再給南塔分配屬性及單元劃分
lsel,s,line,,1,1,1,0
latt,3,25,3,,31,,1
esize,2000
lmesh,all

lsel,s,line,,2,2,1,0
latt,3,25,3,,32,,1
esize,2000
lmesh,all

lsel,s,line,,3,25,2,0
latt,3,25,3,,31,,3
esize,2000
lmesh,all

lsel,s,line,,4,26,2,0
latt,3,25,3,,32,,3
esize,2000
lmesh,all

lsel,s,line,,27,27,1,0
latt,3,25,3,,33,,4
esize,2000
lmesh,all

lsel,s,line,,28,28,1,0
latt,3,25,3,,34,,4
esize,2000
lmesh,all
numstr,elem,4000
SAVE

!ANSYS 在 分網(wǎng)的時候，最后一次操作生成單元之后，相應單元的原來的幾何模型（線、面）并沒有被刪除。還可以顯示。
allsel
eplo

!*****************************************************************

!永和橋的邊界條件

!*****************************************************************
!永和橋的模型邊界條件、成橋荷載、及成橋狀態(tài)靜力分析。尋找邊界分別是有 AXIS里面存儲的梁軸線坐標去尋找的。
! 由于邊界條件對振型和頻率影響十分顯著，所以四個橋墩全部采用彈性邊界條件。然后根據(jù)動力參數(shù)進行修正。
/solu
csys,0

!給北塔北第一個橋墩(橋頭)設置邊界條件（橫向水平約束）
nsel,s,loc,z,axis(1,3)-10,axis(1,3)+10
nsel,r,loc,y,axis(1,2)-10,axis(1,2)+10
nsel,r,loc,x,axis(1,1)+spring-10,axis(1,1)+spring+10
d,all,uy
d,all,ux
d,all,uz

!給北塔北第一個橋墩(橋頭)設置邊界條件（沿橋縱軸向水平約束）
nsel,s,loc,z,axis(1,3)-spring-10,axis(1,3)-spring+10
nsel,r,loc,y,axis(1,2)-10,axis(1,2)+10
nsel,r,loc,x,axis(1,1)-10,axis(1,1)+10
d,all,uy
d,all,ux
d,all,uz
!給北塔北第一個橋墩(橋頭)設置邊界條件（沿橋兩側(cè)豎直方向抗扭轉(zhuǎn)彈簧邊界條件）
nsel,s,loc,z,cabup(1,3)-10,cabup(1,3)+10
nsel,r,loc,y,cabup(1,2)-spring-10,cabup(1,2)-spring+10
nsel,r,loc,x,cabup(1,1)-10,cabup(1,1)+10
d,all,uy
d,all,ux
d,all,uz

nsel,s,loc,z,cabdown(1,3)-10,cabdown(1,3)+10
nsel,r,loc,y,cabdown(1,2)-spring-10,cabdown(1,2)-spring+10
nsel,r,loc,x,cabdown(1,1)-10,cabdown(1,1)+10
d,all,uy
d,all,ux
d,all,uz

!給北塔北第二個橋墩(輔助墩)設置邊界條件（橫向水平約束）
nsel,s,loc,z,axis(5,3)-10,axis(5,3)+10
nsel,r,loc,y,axis(5,2)-10,axis(5,2)+10
nsel,r,loc,x,axis(5,1)+spring-10,axis(5,1)+spring+10
d,all,uy
d,all,ux
d,all,uz

!給北塔北第一個橋墩(輔助墩)設置邊界條件（沿橋兩側(cè)豎直方向抗扭轉(zhuǎn)彈簧邊界條件）
nsel,s,loc,z,cabup(5,3)-10,cabup(5,3)+10
nsel,r,loc,y,cabup(5,2)-spring-10,cabup(5,2)-spring+10
nsel,r,loc,x,cabup(5,1)-10,cabup(5,1)+10
d,all,uy
d,all,ux
d,all,uz

nsel,s,loc,z,cabdown(5,3)-10,cabdown(5,3)+10
nsel,r,loc,y,cabdown(5,2)-spring-10,cabdown(5,2)-spring+10
nsel,r,loc,x,cabdown(5,1)-10,cabdown(5,1)+10
d,all,uy
d,all,ux
d,all,uz

!給北塔設置邊界條件
nsel,s,loc,z,axis(14,3)-10,axis(14,3)+10
nsel,r,loc,y,5000-10,5000+10
d,all,all

!給南塔設置邊界條件
nsel,s,loc,z,axis(38,3)-10,axis(38,3)+10
nsel,r,loc,y,5000-10,5000+10
d,all,all

!給南塔北第二個橋墩(輔助墩)設置邊界條件（橫向水平約束）
nsel,s,loc,z,axis(47,3)-10,axis(47,3)+10
nsel,r,loc,y,axis(47,2)-10,axis(47,2)+10
nsel,r,loc,x,axis(47,1)+spring-10,axis(47,1)+spring+10
d,all,uy
d,all,ux
d,all,uz

!給南塔南第二個橋墩(輔助墩)設置邊界條件（沿橋兩側(cè)豎直方向抗扭轉(zhuǎn)彈簧邊界條件）
nsel,s,loc,z,cabup(47,3)-10,cabup(47,3)+10
nsel,r,loc,y,cabup(47,2)-spring-10,cabup(47,2)-spring+10
nsel,r,loc,x,cabup(47,1)-10,cabup(47,1)+10
d,all,uy
d,all,ux
d,all,uz

nsel,s,loc,z,cabdown(47,3)-10,cabdown(47,3)+10
nsel,r,loc,y,cabdown(47,2)-spring-10,cabdown(47,2)-spring+10
nsel,r,loc,x,cabdown(47,1)-10,cabdown(47,1)+10
d,all,uy
d,all,ux
d,all,uz

!給南塔第一個橋墩(橋頭)設置邊界條件（橫向水平約束）
nsel,s,loc,z,axis(51,3)-10,axis(51,3)+10
nsel,r,loc,y,axis(51,2)-10,axis(51,2)+10
nsel,r,loc,x,axis(51,1)+spring-10,axis(51,1)+spring+10
d,all,uy
d,all,ux
d,all,uz

!給南塔第一個橋墩(橋頭)設置邊界條件（沿橋兩側(cè)豎直方向抗扭轉(zhuǎn)彈簧邊界條件）
nsel,s,loc,z,cabup(51,3)-10,cabup(51,3)+10
nsel,r,loc,y,cabup(51,2)-spring-10,cabup(51,2)-spring+10
nsel,r,loc,x,cabup(51,1)-10,cabup(51,1)+10
d,all,uy
d,all,ux
d,all,uz

nsel,s,loc,z,cabdown(51,3)-10,cabdown(51,3)+10
nsel,r,loc,y,cabdown(51,2)-spring-10,cabdown(51,2)-spring+10
nsel,r,loc,x,cabdown(51,1)-10,cabdown(51,1)+10
d,all,uy
d,all,ux
d,all,uz

!給四個橋墩設置豎向邊界條件（模擬盆式支座）
nsel,s,loc,z,axis(1,3)-10,axis(1,3)+10
nsel,r,loc,y,axis(1,2)-spring-10,axis(1,2)-spring+10
nsel,r,loc,x,axis(1,1)-10,axis(1,1)+10
d,all,uy
d,all,ux
d,all,uz

nsel,s,loc,z,axis(5,3)-10,axis(5,3)+10
nsel,r,loc,y,axis(5,2)-spring-10,axis(5,2)-spring+10
nsel,r,loc,x,axis(5,1)-10,axis(5,1)+10
d,all,uy
d,all,ux
d,all,uz

nsel,s,loc,z,axis(47,3)-10,axis(47,3)+10
nsel,r,loc,y,axis(47,2)-spring-10,axis(47,2)-spring+10
nsel,r,loc,x,axis(47,1)-10,axis(47,1)+10
d,all,uy
d,all,ux
d,all,uz

nsel,s,loc,z,axis(51,3)-10,axis(51,3)+10
nsel,r,loc,y,axis(51,2)-spring-10,axis(51,2)-spring+10
nsel,r,loc,x,axis(51,1)-10,axis(51,1)+10
d,all,uy
d,all,ux
d,all,uz

!給兩個主橋塔處設置垂直橋縱軸線方向水平邊界條件（模擬擋塊約束）
nsel,s,loc,z,axis(14,3)-10,axis(38,3)+10
nsel,r,loc,y,axis(14,2)-10,axis(38,2)+10
nsel,r,loc,x,axis(14,1)+spring-10,axis(38,1)+spring+10
d,all,uy
d,all,ux
d,all,uz

SAVE
allsel
finish

!************************************************************************************************

!動力分析之模態(tài)分析

!************************************************************************************************
allsel

/SOLU
ANTYPE,MODAL
MODOPT,LANB,20
MXPAND,20,,,0
LUMPM,0
PSTRES,on
SOLVE
save
FINISH

/POST1
*GET,OMG1,MODE,1,FREQ
*GET,OMG2,MODE,2,FREQ
*GET,OMG3,MODE,3,FREQ
*GET,OMG4,MODE,4,FREQ
*GET,OMG5,MODE,5,FREQ
*GET,OMG6,MODE,6,FREQ
*GET,OMG7,MODE,7,FREQ
*GET,OMG8,MODE,8,FREQ
*GET,OMG9,MODE,9,FREQ
*GET,OMG10,MODE,10,FREQ

*GET,OMG11,MODE,11,FREQ
*GET,OMG12,MODE,12,FREQ
*GET,OMG13,MODE,13,FREQ
*GET,OMG14,MODE,14,FREQ
*GET,OMG15,MODE,15,FREQ
*GET,OMG16,MODE,16,FREQ
*GET,OMG17,MODE,17,FREQ
*GET,OMG18,MODE,18,FREQ
*GET,OMG19,MODE,19,FREQ
*GET,OMG20,MODE,20,FREQ

*dim,frequence,array,1,20,1

frequence(1,1,1)=omg1
frequence(1,2,1)=omg2
frequence(1,3,1)=omg3
frequence(1,4,1)=omg4
frequence(1,5,1)=omg5
frequence(1,6,1)=omg6
frequence(1,7,1)=omg7
frequence(1,8,1)=omg8
frequence(1,9,1)=omg9
frequence(1,10,1)=omg10
frequence(1,11,1)=omg11
frequence(1,12,1)=omg12
frequence(1,13,1)=omg13
frequence(1,14,1)=omg14
frequence(1,15,1)=omg15
frequence(1,16,1)=omg16
frequence(1,17,1)=omg17
frequence(1,18,1)=omg18
frequence(1,19,1)=omg19
frequence(1,20,1)=omg20

!*cfopen,iteration,dat,'e:\ansyswork\yhdata\output',APPEND,
!*vwrite,NSKIP,NSLOOP(1:1),NSLOOP(1,2),NSLOOP(1,3),NSLOOP(1,4),NSLOOP(2,1),NSLOOP(2,2),NSLOOP(2,3),NSLOOP(2,4),NSLOOP(3,1),NSLOOP(3,2),NSLOOP(3,3),NSLOOP(3,4)
!(F3.0,' ',f8.0,' ',f8.0,' ',f8.0,' ',f8.0,' ',f8.0,' ',f8.0,' ',f8.0,' ',f8.0,' ',f8.0,' ',f8.0,' ',f8.0,' ',f8.0)
!*CFCLOS

*cfopen,optparament,dat,'e:\ansyswork\yhdata\output',APPEND,
*vwrite,NSKIP,optpara(1,1),optpara(1,2),optpara(1,3),optpara(1,4),optpara(1,5),optpara(1,6),optpara(1,7),optpara(1,8),optpara(1,9),optpara(1,10),optpara(1,11),optpara(1,12)
(F3.0,' ',f10.4,' ',f10.4,' ',f10.4,' ',f10.4,' ',f10.4,' ',f10.4,' ',f10.4,' ',f10.4,' ',f10.4,' ',f10.4,' ',f10.4,' ',f10.4)
*CFCLOS

*cfopen,frequence,dat,'e:\ansyswork\yhdata\output',APPEND,
*vwrite,NSKIP,frequence(1,1,1),frequence(1,2,1),frequence(1,3,1),frequence(1,4,1),frequence(1,5,1),frequence(1,6,1),frequence(1,7,1),frequence(1,8,1),frequence(1,9,1),frequence(1,10,1),frequence(1,11,1),frequence(1,12,1),frequence(1,13,1),frequence(1,14,1),frequence(1,15,1),frequence(1,16,1),frequence(1,17,1),frequence(1,18,1)
(F3.0,' ',f10.8,' ',f10.8,' ',f10.8,' ',f10.8,' ',f10.8,' ',f10.8,' ',f10.8,' ',f10.8,' ',f10.8,' ',f10.8,' ',f10.8,' ',f10.8,' ',f10.8,' ',f10.8,' ',f10.8,' ',f10.8,' ',f10.8,' ',f10.8)
*CFCLOS
*SET,object,(OMG1/0.176-1)**2+(OMG2/0.427-1)**2+(OMG3/0.618-1)**2+(OMG6/1.105-1)**2+(OMG11/1.347-1)**2
object=object/5
SAVE

!**************************************************************************************************************************************8

! 提取出模態(tài)位移，保存下來，以供MATLAB 調(diào)用

!****************************************************************************************************************************************

*dim,Beamnode,array,51,1
*dim,towernode,array,4,1
allsel
*do,i,1,51
Beamnode(i,1)=node(axis(i,1),axis(i,2),axis(i,3))
*enddo
towernode(1,1)=node(6250,65000,axis(14,3))
towernode(2,1)=node(-6250,65000,axis(14,3))
towernode(3,1)=node(6250,65000,axis(38,3))
towernode(4,1)=node(-6250,65000,axis(38,3))

! modenum 是選擇存儲的模態(tài)數(shù)
! xmode 是X方向的模態(tài)位移
! ymode 是Y方向的模態(tài)位移
modenum=20
*dim,xmode,array,51,modenum
*dim,ymode,array,51,modenum
*dim,zmode,array,51,modenum
*dim,summode,array,51,modenum
*dim,rotzmode,array,51,modenum
*dim,towermodez,array,4,modenum
*dim,towermodex,array,4,modenum

*do,j,1,modenum
set,1,j
*do,i,1,51
*get,xmode(i,j),node,Beamnode(i,1),u,x
*get,ymode(i,j),node,Beamnode(i,1),u,y
*get,zmode(i,j),node,Beamnode(i,1),u,z
*get,summode(i,j),node,Beamnode(i,1),u,sum
*get,rotzmode(i,j),node,Beamnode(i,1),rot,z
*enddo
*enddo

*do,i,1,51
*do,j,1,modenum
rotzmode(i,j)=rotzmode(i,j)*12000
*enddo
*enddo


*do,j,1,modenum
set,1,j
*do,i,1,4
*get,towermodez(i,j),node,towernode(i,1),u,z
*get,towermodex(i,j),node,towernode(i,1),u,x
*enddo
*enddo
FILENAME='xmode'
*cfopen,%FILENAME%%NSKIP%,dat,'e:\ansyswork\yhdata\output'
*vwrite,sequ,xmode(1,1),xmode(1,2),xmode(1,3),xmode(1,4),xmode(1,5),xmode(1,6),xmode(1,7),xmode(1,8),xmode(1,9),xmode(1,10),xmode(1,11),xmode(1,12),xmode(1,13),xmode(1,14),xmode(1,15),xmode(1,16),xmode(1,17),xmode(1,18)
(F3.0,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12)
*CFCLOS
FILENAME='ymode'
*cfopen,%FILENAME%%NSKIP%,dat,'e:\ansyswork\yhdata\output'
*vwrite,sequ,ymode(1,1),ymode(1,2),ymode(1,3),ymode(1,4),ymode(1,5),ymode(1,6),ymode(1,7),ymode(1,8),ymode(1,9),ymode(1,10),ymode(1,11),ymode(1,12),ymode(1,13),ymode(1,14),ymode(1,15),ymode(1,16),ymode(1,17),ymode(1,18)
(F3.0,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12)
*CFCLOS

FILENAME='zmode'
*cfopen,%FILENAME%%NSKIP%,dat,'e:\ansyswork\yhdata\output'
*vwrite,sequ,zmode(1,1),zmode(1,2),zmode(1,3),zmode(1,4),zmode(1,5),zmode(1,6),zmode(1,7),zmode(1,8),zmode(1,9),zmode(1,10),zmode(1,11),zmode(1,12),zmode(1,13),zmode(1,14),zmode(1,15),zmode(1,16),zmode(1,17),zmode(1,18)
(F3.0,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12)
*CFCLOS
FILENAME='summode'
*cfopen,%FILENAME%%NSKIP%,dat,'e:\ansyswork\yhdata\output'
*vwrite,sequ,summode(1,1),summode(1,2),summode(1,3),summode(1,4),summode(1,5),summode(1,6),summode(1,7),summode(1,8),summode(1,9),summode(1,10),summode(1,11),summode(1,12),summode(1,13),summode(1,14),summode(1,15),summode(1,16),summode(1,17),summode(1,18)
(F3.0,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12)
*CFCLOS
FILENAME='rotzmode'
*cfopen,%FILENAME%%NSKIP%,dat,'e:\ansyswork\yhdata\output'
*vwrite,sequ,rotzmode(1,1),rotzmode(1,2),rotzmode(1,3),rotzmode(1,4),rotzmode(1,5),rotzmode(1,6),rotzmode(1,7),rotzmode(1,8),rotzmode(1,9),rotzmode(1,10),rotzmode(1,11),rotzmode(1,12),rotzmode(1,13),rotzmode(1,14),rotzmode(1,15),rotzmode(1,16),rotzmode(1,17),rotzmode(1,18)
(F3.0,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12)
*CFCLOS
FILENAME='towermodez'
*cfopen,%FILENAME%%NSKIP%,dat,'e:\ansyswork\yhdata\output'
*vwrite,sequ,towermodez(1,1),towermodez(1,2),towermodez(1,3),towermodez(1,4),towermodez(1,5),towermodez(1,6),towermodez(1,7),towermodez(1,8),towermodez(1,9),towermodez(1,10),towermodez(1,11),towermodez(1,12),towermodez(1,13),towermodez(1,14),towermodez(1,15),towermodez(1,16),towermodez(1,17),towermodez(1,18)
(F3.0,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12)
*CFCLOS
FILENAME='towermodex'
*cfopen,%FILENAME%%NSKIP%,dat,'e:\ansyswork\yhdata\output'
*vwrite,sequ,towermodex(1,1),towermodex(1,2),towermodex(1,3),towermodex(1,4),towermodex(1,5),towermodex(1,6),towermodex(1,7),towermodex(1,8),towermodex(1,9),towermodex(1,10),towermodex(1,11),towermodex(1,12),towermodex(1,13),towermodex(1,14),towermodex(1,15),towermodex(1,16),towermodex(1,17),towermodex(1,18)
(F3.0,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12,' ',f16.12)
*CFCLOS