桀outing 发表于 2015-4-20 10:36:04

自编的三维多岩层的强度折减法程序分享,运行成功

n
set log on
set logfile hgtssr.log
rest gravtybalance.sav
;位移清零
ini xdis 0 ydis 0 zdis 0
initial xvelocity 0 yvelocity 0 zvelocity 0
;单元自重应力写入文件'graStress.dat'
def graStressInfo
array buf(6)
pnt=zone_head
loop while pnt # null
    ;应力顺序sxx,syy,szz,sxy,sxz,syz
    buf(1)=string(z_sxx(pnt))
    buf(2)=string(z_syy(pnt))
    buf(3)=string(z_szz(pnt))
    buf(4)=string(z_sxy(pnt))
    buf(5)=string(z_sxz(pnt))
    buf(6)=string(z_syz(pnt))
    status = write(buf,6)
    pnt=z_next(pnt)
end_loop
end
def wrGraStress
status = close
status = open('graStress.dat',1,1)
if status = 0 then
    graStressInfo
    status = close
end_if
end
wrGraStress
;定义调用应力状态的数组,以免重复定义
def defArray
array buu(6)
end


def SSR

;------------------
ait1=0.02         ;精度
k11=1.0          ;上限值
k12=2.0          ;下限值
ks=(k11+k12)/2

loop while (k12-k11)>ait1

;------------------

coh1=5e4/ks
fri1=(atan((tan(20*pi/180))/ks))*180/pi
coh2=7e4/ks
fri2=(atan((tan(28*pi/180))/ks))*180/pi
coh3=7e4/ks
fri3=(atan((tan(28*pi/180))/ks))*180/pi
coh4=9e3/ks
fri4=(atan((tan(18*pi/180))/ks))*180/pi

;------------------
shear1=2.53e8
bulk1=6.60e8
shear2=6.5e8
bulk2=12.6e8
shear3=6.6e8
bulk3=12.6e8
shear4=2.79e8
bulk4=7.28e8
;------------------
   
command

model null       ;挖空
;-----------------------用于读取应力的命令流 --------------------------
;读入单元自重应力
defArray
def graStressInfo
pnt=zone_head
loop while pnt # null
    ;应力顺序sxx,syy,szz,sxy,sxz,syz
    status = read(buf,6)
    z_sxx(pnt)=float(buf(1))
    z_syy(pnt)=float(buf(2))
    z_szz(pnt)=float(buf(3))
    z_sxy(pnt)=float(buf(4))
    z_sxz(pnt)=float(buf(5))
    z_syz(pnt)=float(buf(6))
    pnt=z_next(pnt)
end_loop
end
def reGraStress
status = close
status = open('graStress.dat',0,1)
if status = 0 then
    graStressInfo
    status = close
end_if
end
ini xdisp 0 ydisp 0 zdisp 0
ini xvel 0 yvel 0 zvel 0
model mohr;采用莫尔库仑模型
pro bulk=bulk1 she=shear1 co coh1 fric fri1 ten 0 dilation=10 range group f211 any group f213 any
pro bulk=bulk2 she=shear2 co coh2 fric fri2 ten 8.75e3 dilation=10 range group top any gro 9 any
pro bulk=bulk3 she=shear3 co coh3 fric fri3 ten 8.75e3 dilation=10 range group mid
pro bulk=bulk4 she=shear4 co coh4 fric fri4 ten 3.75e2 dilation=0 range group low
set large

;第1台阶开挖完
m null range gro 9

fix x y z range z 99.9 100.1
fix x range x -0.1 0.1
fix x range x 470.0 472.0
fix y range y -0.1 0.1
fix y range y 265 267
set grav 0 0 -10
set mech ratio 9.8e-5
pl con xdis ou on
solve step 7000
endcommand
;收敛
if mech_ratio<1.0e-4
    k11=ks
    k12=k12
else
    k12=ks
    k11=k11
endif
ks=(k11+k12)/2

endloop
;计算结果的保存
fosfile0='_fos'+'.sav'
command   
save fosfile0
endcommand
end

SSR
pl bl gr
pr ks

桀outing 发表于 2015-4-20 10:37:01

是针对边坡模型的

桀outing 发表于 2015-4-20 10:38:17

断层面的还没进行模拟,只要加上接触面就行了

lijianhao00 发表于 2015-5-13 08:48:30

太棒了 赞一个

zyx0507 发表于 2015-5-13 16:30:34

谢谢分享

周毅 发表于 2015-5-23 17:49:41

学习一下不知道地震荷载怎么施加呢

gyc3017 发表于 2015-5-24 15:16:10

感谢分享

zjt622 发表于 2016-3-10 15:29:27

这个还是不错的

xjtu713 发表于 2016-7-6 18:38:34

好好学习一下!!!

prosper9715 发表于 2016-10-5 10:00:11

fish函数内嵌式定义应该是不行的

企鹅快跑oo 发表于 2016-10-22 13:01:37

很值得学习!谢谢楼主

sjzdh 发表于 2016-10-22 20:24:20

也学习学习

nickjianjyl 发表于 2016-10-22 21:08:02

我也在做这方面的东西,学习学习,请问楼主用的什么建模

414133934 发表于 2016-10-26 13:31:59

感谢楼主分享!!!

爱吃南瓜的然然 发表于 2016-11-9 14:00:46

nickjianjyl 发表于 2016-10-22 21:08
我也在做这方面的东西,学习学习,请问楼主用的什么建模

你知道楼主是用什么建模了么

爱吃南瓜的然然 发表于 2016-11-9 14:02:36

nickjianjyl 发表于 2016-10-22 21:08
我也在做这方面的东西,学习学习,请问楼主用的什么建模

查了一下fish,使用FLAC3D建模的

ahfyclh 发表于 2016-12-2 19:12:59

很值得学习!谢谢楼主

lsj821687340 发表于 2018-8-28 20:16:39

看着不错,学习一下

饿卡酷塔切 发表于 2019-6-6 03:13:58

萨蒂oh地区为hi哦我去hi都饿我去

心无涯 发表于 2019-6-18 21:38:46

为啥我就运行不了啊,你是用的哪个版本软件啊
页: [1] 2
查看完整版本: 自编的三维多岩层的强度折减法程序分享,运行成功