在做桩单元抗滑桩模拟，hist sel pilesel moment 能得出桩身弯矩图吗？我用孙书伟提取梁单元的命令，运行错误，请高手指点一下；还有这段程序能否用在桩单元提取弯矩上？不胜感激！下面附上命令
new
sel beam id=1 begin 0 0 0 end 0 0 10 nseg=20
sel beam prop emod=210e9 nu=0.3 xcarea=1.0 &
xciy=0.0833 xciz=0.0833 xcj=0.167
sel node fix x y z xr yr zr range x=(-0.1,0.1) y=(-0.1,0.1) z=(-0.1,0.1)
sel beam apply zdist 1.0
sel set damp combined
pl sel beam moment my axes
solve
;---------梁单元弯矩输出子程序(基于单元形心)-----
def force_outp
jj=0
loop cid(1,20)
sp=s_find(cid)
if s_type(sp)=1 then
jj=jj+1
;--第1节点相关数据---
np = s_node(sp,1)
xx1 = nd_pos(np,2,1)
yy1 = nd_pos(np,2,2)
zz1 = nd_pos(np,2,3)
ffx1= sb_force(sp,1,1)
ffy1= sb_force(sp,1,2)
ffz1= sb_force(sp,1,3)
mmx1= sb_mom(sp,1,1)
mmy1= sb_mom(sp,1,2)
mmz1= sb_mom(sp,1,3)
;--第2节点相关数据------
np = s_node(sp,2)
xx2 = nd_pos(np,2,1)
yy2 = nd_pos(np,2,2)
zz2 = nd_pos(np,2,3)
ffx2= sb_force(sp,2,1)
ffy2= sb_force(sp,2,2)
ffz2= sb_force(sp,2,3)
mmx2= sb_mom(sp,2,1)
mmy2= sb_mom(sp,2,2)
mmz2= sb_mom(sp,2,3)
; --形心相关数据---
xx=(xx1+xx2)/2.0
yy=(yy1+yy2)/2.0
zz=(zz1+zz2)/2.0
ffx=(ffx1+ffx2)/2.0
ffy=(ffy1+ffy2)/2.0
ffz=(ffz1+ffz2)/2.0
mmx=(mmx1+mmx2)/2.0
mmy=(mmy1+mmy2)/2.0
mmz=(mmz1+mmz2)/2.0
;------------------------
xtable(1,jj)=zz
ytable(1,jj)=mmy
endif
endloop
end
force_outp
set logfile beam_moment.txt
set log on
set pagelength 10000
pri table 1
quit
