ANSYS渡槽槽身動水壓力的施加（1）——矩形渡槽

前言

依據水工抗震規範中關於渡槽動水壓力的部分編一個用於ANSYS渡槽模型動水壓力施加的命令流，是我研究生時一直想要作的一件事，緣由嘛主要是想對比一下規範提供的方法和ANSYS聲學流體單元模擬水體這兩種方法的結果。可是因爲時間關係，畢業之前並無完成，所幸畢業工做以後還有點空閒時間，就把這個命令編了一下。不過我估計幾乎不會有人用這個命令了以及規範的這種方法了，畢竟從理論上來講，用聲學單元模擬獲得的結果確定是比規範的方法精確的，並且操做也很簡單，不涉及編命令。可是無論怎麼樣，仍是放上來充個數吧。

命令依據NB 35047-2015《水電工程水工建築物抗震設計規範》中附錄B 渡槽槽體內動水壓力計算部分編寫，分爲兩個，一個用於矩形渡槽，另外一個用於U型渡槽。

矩形渡槽的動水荷載命令說明

命令運行前的準備工做

• 保證模型的座標系爲：Y方向爲重力方向，X方向爲橫槽向，Z方向爲縱槽向；

• 計算中須要用到ANSYS重啓動技術，由於施加的部分動水荷載的數值須要由渡槽截面槽底中心處的水平或豎直加速度響應值計算獲得，因此命令中加入了導入地震波的命令，須要提早準備好地震波文件accel_x.txt、accel_y.txt、accel_z.txt並記得修改步長DT與總步數NT；
• 需先將槽底面和槽壁面編爲組元：bottom_areas, left_areas,right_areas

• 修改相關參數（紅色）。

pi=3.1415926    !圓周率

H_water=2    !槽內水深

l_water=2    !半槽寬

rou_water=1000   !水體密度

M_water=H_water*l_water*rou_water*2    ! 1m長度的水體質量

node_control=158   !控制點------截面槽底中心處的節點號

 

coord_y_control=ny(node_control)  !獲得控制點的y座標

coord_x_control=nx(node_control)

coord_z_control=nz(node_control)

coord_y_water=coord_y_control+H_water !水面高度座標

 

/prep7

!定義mass21單元

et,99,mass21

!定義彈簧單元

et,98,combin14

!初始實常數

nr=100

 

!讀入地震波

*SET,NT,200      !地震加速度步數

*SET,DT,0.01      !地震加速度步長

*DIM,accel_x,,NT

*DIM,accel_y,,NT

*DIM,accel_z,,NT

*VREAD,accel_x,'accel_x','txt'

(F12.9)

*VREAD,accel_y,'accel_y','txt'

(F12.9)

*VREAD,accel_z,'accel_z','txt'

(F12.9)

 

!***********************************附加質量（開始）**********************************

 !將附加質量的渡槽底面的全部節點編爲一個數組

allsel

CMSEL,S,bottom_areas,AREA

NSLA,S,1

*get,n_bottom,node,,count

*dim,node_bottom,array,n_bottom

*get,nmin,node,,num,min

node_bottom(1)=nmin

*do,i,2,n_bottom

  *get,nnum,node,nmin,nxth

  nmin=nnum

  node_bottom(i)=nmin

*enddo

allsel

 

!將附加質量的渡槽左壁面的全部節點編爲一個數組

CMSEL,S,left_areas,AREA

NSLA,S,1

*get,n_left,node,,count

*dim,node_left,array,n_left

*get,nmin,node,,num,min

node_left(1)=nmin

*do,i,2,n_left

  *get,nnum,node,nmin,nxth

  nmin=nnum

  node_left(i)=nmin

*enddo

allsel

 

!將附加質量的渡槽右壁面的全部節點編爲一個數組

CMSEL,S,right_areas,AREA

NSLA,S,1

*get,n_right,node,,count

*dim,node_right,array,n_right

*get,nmin,node,,num,min

node_right(1)=nmin

*do,i,2,n_right

  *get,nnum,node,nmin,nxth

  nmin=nnum

  node_right(i)=nmin

*enddo

allsel

 

 

!附加渡槽底面節點質量（豎向地震做用衝擊動水壓力）

*do,i,1,n_bottom

nnum=node_bottom(i)

m_wv=0.4*M_water/l_water*ARNODE(nnum)

type,99

nr=nr+1

r,nr,0,0,m_wv

real,nr

e,nnum

*enddo

 

 

!附加渡槽左壁面節點質量（橫向地震做用衝擊動水壓力）

*if,(H_water/l_water),le,1.5,then

*do,i,1,n_left

nnum=node_left(i)

coord_y=ny(nnum)-coord_y_control

m_wh=0.5*M_water/l_water*(coord_y/H_water+0.5*(coord_y/H_water)**2)

m_wh=m_wh*sqrt(3)*tanh(sqrt(3)*l_water/H_water)*ARNODE(nnum)

 

!確保附加的節點位於水面如下

*if,coord_y,le,H_water,then

type,99

nr=nr+1

r,nr,m_wh,0,0

real,nr

e,nnum

*endif

 

*enddo

 

*else

 

*do,i,1,n_left

nnum=node_left(i)

m_wh=0.5*M_water/H_water*ARNODE(nnum)

 

coord_y=ny(nnum)-coord_y_control

!確保附加的節點位於水面如下

*if,coord_y,le,H_water,then

type,99

nr=nr+1

r,nr,m_wh,0,0

real,nr

e,nnum

*endif

 

*enddo

 

*endif

 

!附加渡槽右壁面節點質量（橫向地震做用衝擊動水壓力）

*if,(H_water/l_water),lt,1.5,then

*do,i,1,n_right

nnum=node_right(i)

coord_y=ny(nnum)-coord_y_control

m_wh=0.5*M_water/l_water*(coord_y/H_water+0.5*(coord_y/H_water)**2)

m_wh=m_wh*sqrt(3)*tanh(sqrt(3)*l_water/H_water)*ARNODE(nnum)

 

!確保附加的節點位於水面如下

*if,coord_y,le,H_water,then

type,99

nr=nr+1

r,nr,m_wh,0,0

real,nr

e,nnum

*endif

 

*enddo

 

*else

*do,i,1,n_right

nnum=node_right(i)

m_wh=0.5*M_water/H_water*ARNODE(nnum)

 

coord_y=ny(nnum)-coord_y_control

!確保附加的節點位於水面如下

*if,coord_y,le,H_water,then

type,99

nr=nr+1

r,nr,m_wh,0,0

real,nr

e,nnum

*endif

 

*enddo

*endif

!***********************************附加質量（結束）**********************************

!***********************************施加彈簧（開始）**********************************

g_=9.81

temp=sqrt(5/2)*H_water/l_water

M1=2*rou_water*H_water*l_water*(sqrt(5/2)*l_water/H_water*tanh(temp)/3)

K1=M1*g_/l_water*sqrt(5/2)*tanh(temp)

h1=H_water*(1-(cosh(temp)-2)/temp/sinh(temp))

 

!尋找按高度施加彈簧的節點

coord_y_spring_0=coord_y_control+h1

allsel

CMSEL,S,left_areas,AREA

CMSEL,A,right_areas,AREA

NSLA,S,1

node_find=node(coord_x_control,coord_y_spring_0,coord_z_control)

coord_y_spring_1=ny(node_find)

NSEL,R,loc,Y,coord_y_spring_1-0.001,coord_y_spring_1+0.001

*get,n_spring,node,,count

*dim,node_spring,array,n_spring

*do,i,1,n_spring

  nmin=node(0,0,100000)!選擇距離(0,0,100000)最近的點

  node_spring(i)=nmin

  NSEL,U,node,,nmin!去掉剛選出的節點，從新選擇

*enddo

allsel

*get,node_max,node,,num,max

dist1=0

*do,i,1,n_spring,2

nnum1=node_spring(i)

node1_x=nx(nnum1)

node1_y=ny(nnum1)

node1_z=nz(nnum1)

nnum2=node_spring(i+1)

node2_x=nx(nnum2)

*if,i,eq,n_spring-1,then

dist2=0

*else

nnum3=node_spring(i+2)

node3_z=nz(nnum3)

dist2=node1_z-node3_z

*endif

!建立中間節點

node_max=node_max+1

n,node_max,(node1_x+node2_x)*0.5,node1_y,node1_z

type,99

nr=nr+1

r,nr,M1*0.5*(dist2+dist1),0,0

real,nr

e,node_max

d, node_max,uy

d, node_max,uz

!建立彈簧（模擬對流動水壓力）

type,98

nr=nr+1

r,nr,K1*0.5*(dist2+dist1)

real,nr

e,nnum1,node_max

e,node_max,nnum2

dist1=dist2

*enddo

!***********************************施加彈簧（結束）**********************************

!**********************施加地震加速度和麪荷載形式的動水壓力並計算（開始）******************************

 

!查找槽底槽壁節點對應的渡槽截面中間點（上圖中紅色部分）

allsel

nsel,s,loc,x,coord_x_control-0.001,coord_x_control+0.001

nsel,r,loc,y,coord_y_control-0.001,coord_y_control+0.001

 

*dim,node_bottom_near,array,n_bottom,2

*do,i,1,n_bottom

nnum=node_bottom(i)

nnum_near=nnear(nnum)  !在前面選出的節點中尋找與nnum節點最接近的節點號，並賦值給nnum_near

node_bottom_near(i,1)=nnum_near

*enddo

*dim,node_left_near,array,n_left,2

*do,i,1,n_left

nnum=node_left(i)

nnum_near=nnear(nnum)

node_left_near(i,1)=nnum_near

*enddo

*dim,node_right_near,array,n_right,2

*do,i,1,n_right

nnum=node_right(i)

nnum_near=nnear(nnum)

node_right_near(i,1)=nnum_near

*enddo

 

/solu

ANTYPE,TRANS

TRNOPT,FULL

OUTRES,ALL,ALL

 

allsel

!第一步計算

t=1

time,DT*t

acel,accel_x(t),accel_y(t),accel_z(t)

kbc,0

nsubst,1

rescontrol,DEFINE,last,last !最後時間步重啓

parsav,all,aqueduct_parameter,sav !參數保存

solve

!rescontrol,file_summary

finish

 

*do,t,2,NT

 

!從結果提取控制點加速度

/POST1 

SET,t-1,1

*do,i,1,n_bottom

*get,node_bottom_near(i,2),node,node_bottom_near(i,1),a,x 

*enddo

*do,i,1,n_left

*get,node_left_near(i,2),node,node_left_near(i,1),a,y 

*enddo

*do,i,1,n_right

*get,node_right_near(i,2),node,node_right_near(i,1),a,y

*enddo

parsav,all,aqueduct_parameter,sav !參數保存

 

!再次進入求解模塊

/solu

ANTYPE,,REST,t-1,1,0 

parres,new,aqueduct_parameter,sav !恢復參數

 

!依據渡槽寬高比施加槽底動水壓力荷載（橫向地震做用下槽底的衝擊動水壓力）

*if,(H_water/l_water),lt,1.5,then

*do,i,1,n_bottom

a_wh=node_bottom_near(i,2)

nnum=node_bottom(i)

coord_x_pbh=coord_x_control-nx(nnum)

p_bh=0.5*M_water/l_water*a_wh

p_bh=p_bh*0.5*sqrt(3)*sinh(sqrt(3)*coord_x_pbh/H_water)/cosh(sqrt(3)*l_water/H_water)

F,nnum,FY,p_bh*ARNODE(nnum)

*enddo

*else

*do,i,1,n_bottom

nnum=node_bottom(i)

coord_x_pbh=coord_x_control-nx(nnum)

p_bh=0.5*M_water/H_water/l_water*coord_x_pbh

F,nnum,FY,p_bh*ARNODE(nnum)

*enddo

*endif

 

!左槽壁（豎向地震槽壁動水壓力）

sign_left=(nx(node_left(1))-coord_x_control)/abs(nx(node_left(1))-coord_x_control)

*do,i,1,n_left

a_wv=node_left_near(i,2)

nnum=node_left(i)

coord_y=ny(nnum)-coord_y_control

*if,coord_y,le,H_water,then

p_wv=0.4*M_water/l_water*a_wv*cos(pi*0.5*(coord_y)/H_water)

F,nnum,FX,p_wv*ARNODE(nnum)*sign_left

*endif

*enddo

!右槽壁（豎向地震槽壁動水壓力）

sign_right=(nx(node_right(1))-coord_x_control)/abs(nx(node_right(1))-coord_x_control)

*do,i,1,n_right

a_wv=node_right_near(i,2)

nnum=node_right(i)

coord_y=ny(nnum)-coord_y_control

*if,coord_y,le,H_water,then

p_wv=0.4*M_water/l_water*a_wv*cos(pi*0.5*(coord_y)/H_water)

F,nnum,FX,p_wv*ARNODE(nnum)*sign_right

*endif

*enddo

 

time,DT*t

acel,accel_x(t),accel_y(t),accel_z(t)

nsubst,1

allsel

solve

finish

*enddo

 !**********************施加地震加速度和麪荷載形式的動水壓力並計算（結束）******************************