我編的一個求MC模型單元安全度的FISH程序 mohr-coulomb failure criterion
config zextra 1
def sfdu
p_z=zone_head
loop while p_z # null
c_C=z_prop(p_z,'cohesion')
phi=z_prop(p_z,'friction')
MM1=(z_sxx(p_z)+z_syy(p_z)+z_szz(p_z))
FMAT1=(z_sxx(p_z)-z_syy(p_z))*(z_sxx(p_z)-z_syy(p_z))
FMAT2=(z_syy(p_z)-z_szz(p_z))*(z_syy(p_z)-z_szz(p_z))
FMAT3=(z_szz(p_z)-z_sxx(p_z))*(z_szz(p_z)-z_sxx(p_z))
JIAN1=(z_sxy(p_z))*(z_sxy(p_z))
JIAN2=(z_syz(p_z))*(z_syz(p_z))
JIAN3=(z_sxz(p_z))*(z_sxz(p_z))
MM2=(FMAT1+FMAT2+FMAT3)/6.0+JIAN1+JIAN2+JIAN3
GUANG1=z_sxx(p_z)-MM1/3.0
GUANG2=z_syy(p_z)-MM1/3.0
GUANG3=z_szz(p_z)-MM1/3.0
GUANG4=z_sxy(p_z)*z_syz(p_z)*z_sxz(p_z)
MM3=GUANG1*GUANG2*GUANG3+2.0*GUANG4-GUANG1*JIAN2-GUANG2*JIAN3-GUANG3*JIAN1
GUANG5=sqrt(MM2*MM2*MM2)
kkk=-1*1.5*1.732*MM3/GUANG5
kl=(asin(kkk)/3.0)
k2=sin(kl)
k3=asin(k2)
ss=sqrt(MM2)
SHANG=MM1*sin(phi*pi/180.0)/3.0+(c_C)*cos(phi*pi/180.0)
XIA=cos(kl)+sin(kl)*sin(phi*pi/180.0)/1.732
strengthfactor=(abs(SHANG/XIA))
z_extra(p_z,1)=strengthfactor
p_z=z_next(p_z)
endloop
end
sfdu
pl con zextra 1