Wiki了一下,应变能的定义为应变能密度对体积的积分,应变能密度在各向同性弹性条件下表达式为:
U = 1/2*(sig_x*e_x + sig_y*e_y + sig_z*e_z + tau_xy*gama_xy + tau_yz*gama_yz + tau_zx*gama_zx)
如果用主应力表示为:
U = 1/2(sig1*e1 + sig2*e2 + sig3*e3)
弹性状态下,应力应变一一对应由胡克定律可以把上式改写成:
U = ((1-2v)/6E)*(sig1 + sig2 +sig3)^2 + ((sig1-sig2)^2 + (sig2-sig3)^2 + (sig1-sig3)^2)/12G
E 杨氏模量 v 泊松比 G 剪切模量
上式第一部分为体积应变能,第二部分为偏应变能。
用FISH的话可以设置一个z_extra变量,在计算过程中可以实时显示,代码基于FLAC3D 4.0 或5.0
Def Stain_Energy_elast
E = ....
v = .....
G = E/(2*(1+v))
Sum_U_elast = 0.
p_z = zone_head
loop while zone # null
Uev = ((1-2v)/6E)*(z_sig1(p_z) + z_sig2(p_z) + z_sig3(p_z))^2
Ued = ((z_sig1(p_z) - z_sig2(p_z))^2 + (z_sig2(p_z) - z_sig3(p_z))^2 + (z_sig3(p_z) - z_sig1(p_z))^2)/12/G
z_extra(p_z,1) = (Uev + Ued)*z_volume(p_z)
Sum_U_elast = Sum_U_elast + z_extra(p_z,1)
p_z = z_next(p_z)
endloop
end
@Strain_Energy_elast
hist @Strain_Energy_elast
以上只是考虑到弹性条件下的例子,仅供参考,塑性下需要考虑能量耗散,建议自己查阅文献。
当然如果你用FLAC3D 5.0算得话就不用这么麻烦了,5.0自带了内变量,弹性: z_wevol(ptr), z_weshear(ptr),z_wetot(ptr)
塑性: z_wevol(ptr), z_weshear(ptr),z_wetot(ptr), 直接调用ok啦