用ANSYS画矩形_基于ANSYS的压力自动加载

引言

工程中有一些结构会受到风载荷或水压的作用，如高压线塔或大坝的闸门等。

本篇以一个简单的梁结构为例，演示在ANSYS中如何利用APDL实现梁结构自动加载压力的过程。

Beam188单元

我们经常用Beam188梁单元模拟细长杆结构，如桁架、高压线铁塔、钢架桥等。在ANSYS中自带了如矩形(空心)、圆形(空心)、槽钢、工字型和L型等截面，当然了也可以自定义截面和输入截面常数(面积、转动惯量等)的等效截面，通常情况下用自带的截面形状足以解决大多数工程问题。

在定义Beam188单元的过程中，有一个令人头疼的事，那就是确定截面方向；有些读者可能未处理过该类问题，这里有必要进行简要的解释。

在ANSYS中，Beam188单元的单元坐标系X轴是由节点i指向节点j，Y轴和Z轴分布在截面的两个正交的方向，面载荷(sfbeam、sfe)和单元的应力应变以及轴力弯矩等都是在单元坐标系下定义的。

单元坐标系的Y轴通过截面方向关键点(kp或node)来确定，具体有两种方法：

如果是通过几何线Line离散的方法生成单元，则通过Latt命令中的KB和KE参数确定；
如果是直接建立单元，则在E命令中，第三个节点即方向关键点，如E,1,2,3则单元的ij节点编号分别为1和2，截面方向关键节点k编号为3；

如果不指定方向关键点(kp或node)，则默认单元坐标系的Y轴和总体坐标系的XY平面平行，Z轴则由单元坐标系的XY轴根据右手定则自动确定。只有当截面为圆时，才可不指定方向关键点，否则均需指定方向关键点。

如果指定了方向关键点k(kp或node)，单元坐标系的+X轴方向由i指向j，方向矢量记为Lij，节点i和方向关键点k构成方向矢量Lik，则单元坐标系Y轴为：

即Y轴方向为节点ij和方向关键点k组成平面的法线方向。得到了Y轴的方向之后，Z轴由X轴和Y轴根据右手定则确定。

笔者在ANSYS中尚未找到一个类似*get的命令，让我们可以方便的获取Beam188的单元坐标系，不过我们可根据单元坐标系的定义来获取单元坐标系：

用*get命令得到Beam188单元的节点ij和方向关键点k的三维坐标
计算ij方向(单元坐标系X轴方向)和ik方向的方向矢量Lij和Lik
用*voper命令进行矢量叉乘运算得到单元坐标系Y轴方向矢量
用*voper命令进行矢量叉乘运算得到单元坐标系Z轴方向矢量

以下给出APDL命令流

!获取ijk节点编号
Nodei=nelem(ElemNum,1)
Nodej=nelem(ElemNum,2)
Nodek=nelem(ElemNum,3)
!单元坐标系X轴的方向矢量
ESYSX(1,1)=nx(Nodej)-nx(Nodei)
ESYSX(1,2)=ny(Nodej)-ny(Nodei)
ESYSX(1,3)=nz(Nodej)-nz(Nodei)
!方向矢量Lik
Eik(1,1)=nx(Nodek)-nx(Nodei)
Eik(1,2)=ny(Nodek)-ny(Nodei)
Eik(1,3)=nz(Nodek)-nz(Nodei)
!矢量叉乘运算得到单元坐标系的Y轴和Z轴
*voper,ESYSY,Eik,cross,ESYSX
*voper,ESYSZ,ESYSX,cross,ESYSY

注意：在ANSYS中进行矢量的叉乘运算用*voper，关键字cross，矢量必须用一维行数组存放，即1行3列的数组，而不能用3行1列的数组

Beam188单元的面载荷通过sfbeam命令施加，面载荷种类分别为上图中的12345，为了叙述方便，假设力、长度和压力单位分别为N、m和Pa，当然也可根据实际情况定为N、mm和Mpa。

1表示单元X轴方向的面1受到单位长度的力作用，即面载荷单位为N/m,正方向为-Z方向，如截面为矩形，则面1表示梁的上顶面；
2表示单元X轴方向的面2受到单位长度的力作用，即面载荷单位为N/m,正方向为-Y方向，如截面为矩形，则面2表示梁的右侧面；
3表示单元X轴方向受到单位长度的力作用，即面载荷单位为N/m，正方向为+X方向，面3表示梁的长度方向；
4表示单元i节点所在的截面受到的集中力作用，即单位为N，正方向为+X方向；
5表示单元j节点所在的截面受到的集中力作用，即单位为N，正方向为-X方向；

注意：前三种面载荷的单位为N/m，后两种载荷的单位为N，当外载荷为风压或者水压时，我们一般得到的输入数据是单位为Pa的压强，因此，当施加前三种面载荷时，应把风压或水压乘以梁截面的宽度或高度，即Pa*m=N/m；当施加后两种面载荷时，应该把面载荷值乘以横截面面积，即Pa*m2=N。

典型命令

SFBEAM, Elem, LKEY, Lab, VALI, VALJ, VAL2I, VAL2J, IOFFST, JOFFST, LENRAT

——在Beam单元上加面载荷

注意：IOFFST和JOFFST选项只能在Beam188单元的Keyopt(3)=3，即形函数必须为三次时才可用。

算例

我们在做有限元分析时，会在Beam188梁单元建模的结构上施加如风压或水压等面载荷，需要用sfbeam命令进行加载，不同的面因为尺寸不一致，所以受到的面力也是不一样的，因此需要根据外载荷方向和单元坐标系的方向确定具体哪一个面受载，具体来说就是要把外载荷沿单元坐标系分解，根据单元坐标系确定受载面。

下面以Beam188梁单元建模的简单框架结构为例，结构受到水压作用，水面位于Y=0位置，演示如何实现水压的自动加载。

以下为自动加载水压的命令流，限于篇幅，建模部分的命令流就不在此列出，有兴趣的读者可私信交流获取，值得注意的是，把外载荷分别和单元坐标系XYZ轴的方向矢量分别做点乘，然后除以XYZ轴方向矢量的模，即可得到外载荷在单元坐标系下的分量，最后再加载即可，就变成了加风载荷的命令流。

finish
/clear
/psf,pres,norm,2,,on
/psymb,esys,on*set,DensWater,1000
*set,GravityVal,10
*dim,PDir,array,1,3,1
PDir(1,1)=-1
PDir(1,2)=0
PDir(1,3)=-1
*dim,ECSYSX,array,1,3,1
*dim,ECSYSY,array,1,3,1
*dim,ECSYSZ,array,1,3,1
*dim,Eik,array,1,3,1 /prep7
!建模开始
! ......
!建模结束
*get,ElemCount,elem,0,count
*get,ElemMin,elem,0,num,min
*get,ElemMax,elem,0,num,max
*dim,ElemNum,array,ElemCount,1,1
*dim,VDot,array,1,1,1
!获取单元编号
j=0
*do,i,ElemMin,ElemMax,1 *if,esel(i),eq,1,thenj=j+1ElemNum(j)=i*endif
*enddo!循环单元加载
*do,i,1,ELemCount,1!单元节点编号          Nodei=nelem(ElemNum(i),1) Nodej=nelem(ElemNum(i),2) Nodek=nelem(ElemNum(i),3) !单元坐标系X轴方向矢量ECSYSX(1,1)=nx(Nodej)-nx(Nodei)   ECSYSX(1,2)=ny(Nodej)-ny(Nodei)    ECSYSX(1,3)=nz(Nodej)-nz(Nodei)    !i节点到方向关键点k的方向矢量Eik(1,1)=nx(Nodek)-nx(Nodei)   Eik(1,2)=ny(Nodek)-ny(Nodei)   Eik(1,3)=nz(Nodek)-nz(Nodei)!矢量叉乘运算得到单元坐标系的YZ轴*voper,ESYSY,Eik,cross,ESYSX *voper,ESYSZ,ECSYSX,cross,ESYSY!获取梁截面尺寸*get,SecNum,elem,ElemNum(i),attr,secn *get,ElemSecY,secp,SecNum,data,1   *get,ElemSecZ,secp,SecNum,data,2    !水压方向矢量和单元坐标系Z轴点乘*voper,VDot,PDir,dot,ESYSZ *if,VDot(1,1),GE,1.0E-8,then
Pi= -DensWater*GravityVal*abs(ny(Nodei))*ElemSecY
Pj=-DensWater*GravityVal*abs(ny(Nodej))*ElemSecY
sfbeam,ElemNum(i),1,press,Pi,Pj*endif*if,VDot(1,1),LE,1.0E-8,then
Pi= DensWater*GravityVal*abs(ny(Nodei))*ElemSecY
Pj=DensWater*GravityVal*abs(ny(Nodej))*ElemSecY
sfbeam,ElemNum(i),1,press,Pi,Pj*endif
!水压方向矢量和单元坐标系Y轴点乘*voper,VDot,PDir,dot,ESYSY *if,VDot(1,1),GE,1.0E-8,then
Pi= -DensWater*GravityVal*abs(ny(Nodei))*ElemSecZ
Pj=-DensWater*GravityVal*abs(ny(Nodej))*ElemSecZ
sfbeam,ElemNum(i),2,press,Pi,Pj*endif*if,VDot(1,1),LE,1.0E-8,then
Pi= DensWater*GravityVal*abs(ny(Nodei))*ElemSecZ
Pj=DensWater*GravityVal*abs(ny(Nodej))*ElemSecZ
sfbeam,ElemNum(i),2,press,Pi,Pj*endif
*enddo
finish

-完-

公众号同步更新：数联科技工作室

公众号内有更多内容！