PFC2D5.0_双轴试验

R_S^^

已于 2024-01-18 16:02:42 修改

阅读量1.5k

点赞数 19

分类专栏： PFC 文章标签：经验分享

于 2024-01-18 15:32:13 首次发布

本文链接：https://blog.csdn.net/qq_56811637/article/details/135669587

版权

PFC 专栏收录该内容

6 篇文章

订阅专栏

本文详细介绍了使用PFC软件进行土体模拟实验，包括等向压缩和双向剪切过程中的代码片段，展示了如何通过伺服机制控制墙体运动并调整力，以达到预定的围压。作者也鼓励读者就相关土木工程或岩土问题进行讨论。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

一、效果展示
Specimen generation
Isotropic compression（本文施加围压500kPa）
在这里插入图片描述

Biaxial shearing
在这里插入图片描述

二、代码展示

new
set random 10000
;-------------------------------Specimen generation--------------------------------
[L=0.035]
[H=0.035]
domain extent [-0.8*L] [0.8*L] [-0.8*H] [0.8*H] condition stop 
;;;; Create domain
wall generate name 'box' box [-0.5*L] [0.5*L] [-0.5*H] [0.5*H] expand 1.5 
;;;; Generate wall
ball distribute porosity 0.15 radius 0.0003 0.00045 ...
     box [-0.5*L] [0.5*L] [-0.5*H] [0.5*H]

ball attribute density 1000 damp 0.7
;;;; Set ball density and damping ratio
cmat default model linear prop kn 5e8 ks 5e8 fric 0.3 type 'ball-ball'
cmat default model linear prop kn 5e8 ks 5e8 fric 0.3 type 'ball-facet'
;pause key
;;;; Set contact model and assign model parameters
cycle 5000 calm 10
;pause key
cmat default model linear prop kn 5e8 ks 5e8 fric 0.3 type 'ball-ball'
cmat default model linear prop kn 5e7 ks 5e7 fric 0.3 type 'ball-facet'
cycle 1000
;wall list
pause key
;;-------------------------------Isotropic compression---------------------------
;;define servo_walls
;;      wall.servo.force.x(wall.find('boxRight')) =   -txx*wly
;;      wall.servo.force.x(wall.find('boxLeft'))	=    txx*wly
;;      wall.servo.force.y(wall.find('boxBack'))	=   -tyy*wlx
;;      wall.servo.force.y(wall.find('boxFront'))	=    tyy*wlx
;;end
;;set fish callback  9.0 @servo_walls
define stress_disp
wlx=wall.pos.x(wall.find(2))-wall.pos.x(wall.find(4))
wly=wall.pos.y(wall.find(3))-wall.pos.y(wall.find(1))
wsyy=(wall.force.contact.y(wall.find(3))-wall.force.contact.y(wall.find(1)))/wlx*0.5
wsxx=(wall.force.contact.x(wall.find(2))-wall.force.contact.x(wall.find(4)))/wly*0.5
end
@stress_disp

[txx=0.5*1.0e6]
[tyy=0.5*1.0e6]
wall servo activate on yforce [-tyy*wlx] vmax 1.00 range set name 'boxTop'
wall servo activate on yforce [tyy*wlx] vmax 1.00 range set name 'boxBottom'
wall servo activate on xforce [-txx*wly] vmax 1.00 range set name 'boxRight'
wall servo activate on xforce [txx*wly] vmax 1.00 range set name 'boxLeft'
define servo_walls
      stress_disp
      wall.servo.force.x(wall.find(2)) 	=   -txx*wly
      wall.servo.force.x(wall.find(4))	=    txx*wly
      wall.servo.force.y(wall.find(3))	=   -tyy*wlx
      wall.servo.force.y(wall.find(1))	=    tyy*wlx
end
set fish callback  9.0 @servo_walls

hist id 1 @wsxx
hist id 2 @wsyy
plot create plot 'wsxx vs step'
plot add hist 1 2 vs step

define iso_compression

tol=0.005
loop while 1#0
  if math.abs((wsyy - tyy)/tyy) < tol
    if math.abs((wsxx-txx)/txx)<tol
       exit
    endif
  endif
   command
    cycle 100
   endcommand
endloop

end
@iso_compression
pause key
;-------------------------------Biaxial shearing------------------------------------------
[rate = 0.2]
[wly0=wly]
[weyy=0]
[wexx=0]
[vol=0]
[wlx0=wlx]
[devi=wsyy-txx]
hist id 3 @weyy
hist id 4 @devi
hist id 5 @vol
plot create plot 'Deviatoric stress vs Axial strain'
plot add hist 4 vs 3
plot create plot 'Volumetric strain vs Axial strain'
plot add hist 5 vs 3

ball attribute displacement multiply 0.0
wall attribute displacement multiply 0.0
;pause key
wall servo activate off range set name 'boxTop' set name 'boxBottom' union
;;;;Close the servo mechanism of the top and bottom walls 
wall attribute yvelocity [-rate*wly0] range set name 'boxTop'
wall attribute yvelocity [ rate*wly0] range set name 'boxBottom'
;;;;Assign a constant loading velocity to the Top and Bottom wall
[stop_me = 0]
[target = 0.075]
define stop_me
    devi=wsyy-txx
    ;weyy=math.abs(wall.disp.y(wall.find(3))-wall.disp.y(wall.find(1)))/wly0
     weyy=(wly0-wall.pos.y(wall.find(3))+wall.pos.y(wall.find(1)))/wly0
     wexx=(wlx0-wall.pos.x(wall.find(2))+wall.pos.x(wall.find(4)))/wlx0
     vol=-weyy-wexx
  if weyy > target then
    stop_me = 1
  endif
end
solve fishhalt @stop_me

三、关于伺服的讨论
开启墙体伺服机制：

[txx=0.5*1.0e6]
[tyy=0.5*1.0e6]
wall servo activate on yforce [-tyy*wlx] vmax 1.00 range set name 'boxTop'
wall servo activate on yforce [tyy*wlx] vmax 1.00 range set name 'boxBottom'
wall servo activate on xforce [-txx*wly] vmax 1.00 range set name 'boxRight'
wall servo activate on xforce [txx*wly] vmax 1.00 range set name 'boxLeft'

定义及调用伺服函数：

define servo_walls
      stress_disp
      wall.servo.force.x(wall.find(2)) 	=   -txx*wly
      wall.servo.force.x(wall.find(4))	=    txx*wly
      wall.servo.force.y(wall.find(3))	=   -tyy*wlx
      wall.servo.force.y(wall.find(1))	=    tyy*wlx
end
set fish callback  9.0 @servo_walls

伺服收敛条件：

define iso_compression

tol=0.005
loop while 1#0
  if math.abs((wsyy - tyy)/tyy) < tol
    if math.abs((wsxx-txx)/txx)<tol
       exit
    endif
  endif
   command
    cycle 100
   endcommand
endloop

end

在加载过程中，墙体向土体方向移动，wlx、wly在一直发生变化，为达到指定的围压，wall servo命令指定的力必须随着wlx、wly的变化而调整，调用伺服函数和callback命令就可以实现这个目的。iso_compression函数的调用是判断围压是否达到500kPa（误差0.005），若达到，则结束。

有任何PFC问题（土木方向，最好岩土），欢迎与我讨论，一起解决问题。