工程实例
1、FLAC3D 锚杆建模中出现的问题
请教各位大虾：
我用FLAC3d建模计算锚杆的应力和位移，发现几个问题不知如何处理？
1. 计算出来的位移值偏大，达到10cm，而施加
荷载还未达到设计值。
2. 锚杆周边的位移和应力云图好像不太正常。 new Title
gen zone brick p0 0 0 0 p1 6 0 0 p2 0 0 -12 p3 0 6 0 size 15 30 15 ratio 1.0 1.0 1.0
gen zone reflect normal 0 -1 0 origin 0 0 0 gen zone reflect normal -1 0 0 origin 0 0 0 model mohr
prop bulk 9.52e6 shear 4.65e6 coh 0 fric 25 ini dens 2000
fix z range z -11.9 -12.1 fix y range y 5.9 6.1 fix y range y -6.1 -5.9 fix x range x 5.9 6.1 fix x range x -6.1 -5.9 set grav 0 0 -10 set large hist unbal
sel cable id 1 begin 0.037 0 0.9 end 0.037 0 -3.6 nseg 15
sel cable id 2 begin -0.037 0 0.9 end -0.037 0 -3.6 nseg 15
sel cable prop xcarea 0.0167 emod 200e9 ytens 335e6 gr_k 16.7e6 gr_per 0.459 gr_coh 200
sel cable id 1 pre 5e4 range x 0.0369 0.0371 y -0.001 0.001 z 0.69 0.71
sel cable id 2 pre 5e4 range x -0.0371 -0.0369 y -0.001 0.001 z 0.69 0.71
plot set rotation 30 0 40 plot set center auto
plot add contour disp outline on plot add sketch plot add disp
plot add sel cable force red yellow plot add sel cable grout stress blue plot add sel cable grout slip green plot show step 2000
2、锚杆受拉超限
锚杆能够承受拉力，而不能承受超过规定的剪应力，岩体的变形超过了锚杆抗剪范围，会出现上述警告。 锚杆的端头是否经过预应力锚杆处理的，关于预应力锚杆的模拟，你可以到仿真论坛 lakerwater的帖子看看。 ;计算中锚杆的长度发生了变化!!!!!!!!!!!!!
;如果不施加初锚力，锚杆的长度不会变，但是一旦施加了初锚力，锚杆的长度
;就会发生变化，能直接从图中看出，并且越靠近隧道口（y=0附近），锚杆长度
;变化越明显。隧道尾（y=20附近）处锚杆的长度变化不明显。
;同时，当初锚力达到6t时，锚杆变化异常，出现了很多警告。如下：
;+++ WARNING: axial dir. of cables and piles or normal dir. of liners
;and geogrids at shared node out of alignment by more than 30 degrees.
;这是怎么回事？请大侠赐教。。。。。。。
rest Model.sav ;施加2t的初锚力 set log on
set logfile 锚杆2t.log ;材料模式：摩尔-库伦 m m
;========围岩参数========= ;加载层(细砂岩)
pro b 10.9e8 s 8.7e8 f 40 c .93e5 t .93e4 ran z 47.12 70 ini d 2500 ran z 47.12 70 ;-----------------------
;老顶(中粒砂岩)
pro b 12.1e9 s 10.8e9 f 44 c 1.18e6 t 1.18e5 ran z 36.91 47.12
ini d 2500 ran z 36.91 47.12 ;----------------------- ;直接顶(细砂岩)
pro b 10.9e8 s 8.7e8 f 40 c .93e5 t .93e4 ran z 33.32 36.91 ini d 2500 ran z 33.32 36.91 ;----------------------- ;伪顶(碳质泥岩)
pro b 7.1e8 s 4.7e8 f 33 c .53e5 t .53e4 ran z 32.8 33.32 ini d 2500 ran z 32.8 33.32 ;----------------------- ;煤层(二号煤层)
pro b 7.1e8 s 4.7e8 f 33 c .53e5 t .53e4 ran z 30 32.8 ini d 2300 ran z 30 32.8
21
;-----------------------
;直接底(泥岩)
pro b 9.3e8 s 7.1e8 f 37 c .77e5 t .77e4 ran z 25 30 ini d 2500 ran z 25 30 ;-----------------------
;老底(细砂岩)
pro b 11.8e9 s 10.1e9 f 43 c 1.09e6 t 1.09e5 ran z 0 25 ini d 2500 ran z 0 25
;=============================================
;固定边界
fix x ran x -0.1 0.1 fix x ran x 74.9 75.1 fix y ran y -0.1 0.1 fix y ran y 19.9 20.1 fix z ran z -0.1 0.1
;=============================================
;原岩应力场
ini sxx -25.75e6 grad 0 0 2.5e4 ran z 0 70 ini syy -25.75e6 grad 0 0 2.5e4 ran z 0 70 ini szz -25.75e6 grad 0 0 2.5e4 ran z 0 70 app szz -24e6 ran z 70
;======================================== ;计算设置
set gravity 0 0 -10 set large
;====================================== ;显示垂直应力场分布 pl bcon szz out on
;====================================== ;原岩应力场平衡计算 solve
;====================================== ;调整生成原岩应力场后的节点速度和位移为0 ini xvel 0.0 yvel 0.0 zvel 0.0 ini xd 0.0 yd 0.0 zd 0.0
;====================================== ;开挖巷道
m n ran gro 10 an gr 11 an gr 12 an gr 13 an
;=====定义锚杆，长度2m=================== def inbolts
xtable(1,1)=39.7524 xtable(1,2)=40
xtable(1,3)=36.9437 xtable(1,4)=38.0563 xtable(1,5)=39.0587 xtable(1,6)=35.9413
xtable(1,7)=35
xtable(1,8)=35.2476 xtable(2,1)=41.5544 xtable(2,2)=42
xtable(2,3)=36.4987 xtable(2,4)=38.5013 xtable(2,5)=40.3057 xtable(2,6)=34.6943 xtable(2,7)=33
xtable(2,8)=33.4456
ytable(1,1)=28.5847 ytable(1,2)=27.5 ytable(1,3)=29.9373 ytable(1,4)=29.9373 ytable(1,5)=29.4546 ytable(1,6)=29.4546 ytable(1,7)=27.5 ytable(1,8)=28.5847 ytable(2,1)=29.4525 ytable(2,2)=27.5 ytable(2,3)=31.8872 ytable(2,4)=31.8872 ytable(2,5)=31.0182 ytable(2,6)=31.0182 ytable(2,7)=27.5 ytable(2,8)=29.4525
nn=0 ;锚杆单元序列号 ;隧道锚杆 loop ii (1,19)
loop jj (1,8) nn=nn+1
bx=xtable(1,jj) bz=ytable(1,jj) ex=xtable(2,jj) ez=ytable(2,jj) command
sel cable id nn beg bx ii bz end ex ii ez nseg 4
sel cable pre 60000 prop emod 45e9 xcarea 3.8e-4 gr_per 1 yten 2.5e6 gr_k 17.5e7 gr_c 2e6 endcommand endloop endloop end inbolts
;======================================== ;显示地层移动分布 pl con zdisp out on
22
pl ad sel geo scale 0.005
pl set center 36.7 7.814 28.61 pl set rotation 20 0 20 pl set mag 9.31 pl set dist 217.4 print sel ca le
;======================================== ;巷道开挖平衡计算 solve
print sel ca le sav 4-6.sav
3、FLAC3D计算完毕关机
4、隧道模拟
参数：;问题描述 考虑圣维南原理,取周围岩土的尺寸维隧道尺寸的5~6倍,此处取为15m r;初衬 C50管片外径3.0m，内径2.7m，C50,E＝34.5GPa，V＝0.167,.g;二衬 外径2.7m，内径2.3m,C30,E＝30GPa，V＝0.167（二衬先不考虑）& ;粘质粉土，4.8m，φ＝26°，重度＝19，c＝11KPa，变形模量＝10MPa
;中细纱，1.7m，φ＝28°，重度＝20，c＝0KPa，变形模量＝13Mpa) ;问题描述 考虑圣维南原理,取周围岩土的尺寸维隧道尺寸的5~6倍,此处取为15m
;初衬 C50管片外径3.0m，内径2.7m，C50,E＝34.5GPa，V＝0.167,
;二衬 外径2.7m，内径2.3m,C30,E＝30GPa，V＝0.167
;粘质粉土，4.8m，φ＝26°，重度＝19，c＝11KPa，变形模量＝10MPa+ |& c. g% T P w
;中细纱，1.7m，φ＝28°，重度＝20，c＝0KPa，变形模量＝13Mpa
;-------------------------------前处理 title '沉降计算'
plot set title text 'chenjiang'8 d0 v0 b* o% t plot set rotation 20 0 30 5 e/ \\: b' o. u6 |: Y plot set center 0 0 5 plot set dist 300 plot set mag 0.8
plot add surface green
plot add axes0 @0 G( b c- Q7 D, M, s2 f5 _ plot show
;建立模型---------------------------------4 d gen zone radcyl p0 0 0 0 p1 4.7 0 0 p2 0 50 0 p3 0 0 4.7 & size 10 30 10 4 dim 3 3 3 3 group tunnel gen zone brick p0 4.7 0 0 p1 18 0 0 p2 4.7 50 0 p3 4.7 0 4.7 &size 10 30 5 group brick1! gen zone cshell p0 0 0 0 p1 3 0 0 p2 0 50 0 p3 0 0 3 & dim 2.7 2.7 2.7 2.7 size 1 30 10 6 group segment gen zone cylinder p0 0 0 0 p1 2.7 0 0 p2 0 50 0 p3 0 0 2.7 &*size 6 30 10 group cy, s6 X' K; Z' k/ r( E
gen zone reflect ori 0 0 0 norm 0 0 -1 range z 0 4.73 gen zone brick p0 0 0 4.7 p1 4.7 0 4.7 p2 0 50 4.7 p3 0 0 9.5 &)size 5 30 5 group brick2) X gen zone brick p0 4.7 0 4.7 p1 18 0 4.7 p2 4.7 50 4.7 p3 4.7 0 9.5 &#size 10 30 5 group brick3. gen zone brick p0 0 0 -18 p1 4.7 0 -18 p2 0 50 -18 p3 0 0 -4.7 &$size 5 30 15 group brick4gen zone brick p0 4.7 0 -18 p1 18 0 -18 p2 4.7 50 -18 p3 4.7 0 -4.7 & size 10 30 15 group brick5# gen zon reflect norm 1 0 0 orig 0 0 0
;-赋土参数-------------------------------------------4 model mohr ;E=13Mpa,v=0.3
prop density 2000 bulk 13e6 shear 5eric 28 coh 1e10 tension 1e10 range z -18 4.7
;E=10Mpa,v=0.255 E, s% d) u ?# J( _+ B
prop density 1900 bulk 10e6 shear 4eric 26 coh 1e10 tension 1e10 range z 4.7 9.5$ X ;初始地应力--------------------------------------! ini xdisp 0 ydisp 0 zdisp 0
ini szz -1.85e5 grad 0 0 1.9e4 range z 4.7 9.5 ;19×4.8+20×4.7＝185.2KPa) |' M8 F' c+ E+ d ini sxx -6.67e4 grad 0 0 6.84e3 range z 4.7 9.5 ;侧压力系数取为0.36 （66.67KPa）
ini syy -6.67e4 grad 0 0 6.84e3 range z 4.7 9.5* k) ~8 `0 q% g ini szz -1.85e5 grad 0 0 2.0e4 range z -18 4.7
23
ini sxx -6.67e4 grad 0 0 7.2e3 range z -18 4.7 ini syy -6.67e4 grad 0 0 7.2e3 range z -18 4.7 ;边界条件-------------------------------------
apply nstress -1e4 range z 9.4 9.6* B( h+ C. o! m5 @! D
fix x range x -18.1 -17.9
fix x range x 17.9 18.1/ s1 W( J0 Y& U2 I% c& j fix y range y 49.9 50.1
fix y range y -0.1 0.1# F$ O8 X( Z. V- e1 Y0 n fix z range z -18.1 -17.9 set grav 0 0 -10 small
;solve. u7 |. @$ O9 ^7 r' K0 H' d
;save small.sav0 \\$ R+ S3 l, W$ u8 f\ - |% ]' n# @; `; Q! y ;开挖 求解
model null range group cy/ b r;---------------------------! ini xdisp 0 ydisp 0 zdisp 0 ;---------------------------
model elas range group segment
prop density 2500 bulk 34.5e9 shear 14.78e9 range group segment
;E=13Mpa,v=0.3* n$ {5 D9 j$ q m3 O7 ?+ u$ J: w3 prop density 2000 bulk 13e6 shear 5eric 28 coh 0 range group brick1
prop density 2000 bulk 13e6 shear 5eric 28 coh 0 range group tunnel/ b$ z' m/ V: d6 M& A' t! S
prop density 2000 bulk 13e6 shear 5eric 28 coh 0 range group brick4
prop density 2000 bulk 13e6 shear 5eric 28 coh 0 range group brick5
;E=10Mpa,v=0.251 j; W9 @& N7 o ~+ v3 G
prop density 2000 bulk 13e6 shear 5eric 28 coh 0 range z 4.7 9.5\hist unbal
hist gp xdis 3 0 0& o, B: N! @\ hist gp xdis -3 0 0 hist gp zdis 0 0 9.5
hist gp zdis 0 0 7.51 y1 N) h: N( `' u hist gp zdis 0 0 5.5. t\ hist gp zdis 0 0 4.5 hist gp zdis 0 0 3
hist gp zdis 0 0 -3* Z8 t, O' u5 Z1 @6 z set large- Q1 ?8 p. w7 E) L\ ;solve ;save n.sav
5、显示最大位移点的坐标参考
def get_gp_maxdisp3 R2 R) n) M( Q# f! D$ z
gp0_disp = gp_xdisp(gp_head)*gp_xdisp(gp_head) gp0_disp = gp0_disp + gp_ydisp(gp_head)*gp_ydisp(gp_head) gp0_disp = gp0_disp + gp_zdisp(gp_head)*gp_zdisp(gp_head) gp0_disp = sqrt(gp0_disp) # S/ N+ h( @# C( E& }0 p- O% p_gp=gp_head: ^* E @6 e6 Z$ r ;找最大值* H/ ~4 `. m; ?/ E
loop while p_gp # null8 ]( v9 N* \\9 K) L l; X2 M gp_disp = gp_xdisp(p_gp)*gp_xdisp(p_gp)
gp_disp = gp_disp + gp_ydisp(p_gp)*gp_ydisp(p_gp) gp_disp = gp_disp + gp_zdisp(p_gp)*gp_zdisp(p_gp) gp_disp = sqrt(gp_disp) : U* q! w+ w& k5 ?/ g if gp_disp>gp0_disp gp0_disp=gp_disp endif ' Y' r( ` p_gp = gp_next(p_gp)% Z3 i4 b4 t$ H endloop ;找最大值的坐标 p_gp=gp_head ss=0+ ]1 loop while p_gp # null
gp_disp = gp_xdisp(p_gp)*gp_xdisp(p_gp)) I7 D Z0 z gp_disp = gp_disp + gp_ydisp(p_gp)*gp_ydisp(p_gp) gp_disp = gp_disp + gp_zdisp(p_gp)*gp_zdisp(p_gp) gp_disp = sqrt(gp_disp)
if gp_disp=gp0_disp. u5 [( r$ H/ i' m ss=ss+1
xtable(1,ss)=gp_xpos(p_gp) ytable(1,ss)=gp_ypos(p_gp) endif
p_gp = gp_next(p_gp) endloop end
get_gp_maxdisp
6、固结问题输出某一节点随时间的孔压变化
方法其实很简单，如果不需要导出数据，直接hist gp pp就可以，如果要导出数据的话就需要使用fish。例如你要记录每隔1hour即3600s后某一节点的孔压，首先你要用fish来控制solve的时间间隔，然后把每个小时后的孔压记录到一个table中，最后plot table就可以了，要把数据导出时，现log on，然后print table(),然后log off即可，具体操作如下,假设求解10小时，所记录节点坐标为0，0，0：
def pplog
sage=3600 ;渗流时间 % p, Y$ C& F4 `- [9 i. A\ pnt=gp_near(0,0,0) n=1
loop while sage <= 36000 ;计算时间以10小时
24
