单向流固耦合——模拟弯曲河床砂石运移

lobby 2022年7月28日浏览：2304 评论：2 收藏：3

0 引言

要完成这个工况的模拟，其实需要解决两个难点，第一个是如何生成河床。这里假设河床截面是圆形的，那其实是需要一个pipe形状的三维体，但是这个在PFC中是没有提供的。第二个难点就是河流流场的生成了。

1 生成河床

这里我们使用的是以直代曲的办法，可以用一个个圆柱拼接成河床的额样子。createwall 函数中就是使用循环生成河床的墙体，R_he 为河床截面中心的弧长半径，jing_he 为河床截面的圆形半径，he_kua 为河床跨过的圆心角。wallsplit是我们河床墙体的分段，也不是说分的越细结果越好的。

ball_range 我们定义了河床砂石的生成区域，也是一个圆柱形区域，用的是generate生成后自重沉降的方式。

new domain extent 0 20 -20 20 -4 4

[wallsplit=50][R_he=16][jing_he=2]

[he_kua=math.pi*1.0/3.0]def createwall    jiaodu_split=he_kua/float(wallsplit)    loop n(1,wallsplit)        jiaodu1=-he_kua*0.5+jiaodu_split*(n-1)        jiaodu2=-he_kua*0.5+jiaodu_split*(n)        x_pos1=math.cos(jiaodu1)*R_he        y_pos1=math.sin(jiaodu1)*R_he        x_pos2=math.cos(jiaodu2)*R_he        y_pos2=math.sin(jiaodu2)*R_he        vec=math.unit(vector(x_pos2-x_pos1,y_pos2-y_pos1,0))        height=(R_he+jing_he)*jiaodu_split        command            wall generate cylinder base @x_pos1 @y_pos1 0 axis @vec height @height radius @jing_he ...                cap false false resolution  0.5             endcommand    endloopend@createwall wall delete range z 0 20 extent

def ball_range    jiaodu1=-he_kua*0.5+jiaodu_split*6    jiaodu2=-he_kua*0.5+jiaodu_split*7    x_pos1=math.cos(jiaodu1)*R_he    y_pos1=math.sin(jiaodu1)*R_he    x_pos2=math.cos(jiaodu2)*R_he    y_pos2=math.sin(jiaodu2)*R_he    vec=math.unit(vector(x_pos2-x_pos1,y_pos2-y_pos1,0))    height=(R_he+jing_he)*jiaodu_split*8    command        geometry set ballrange        geometry generate cylinder base @x_pos1 @y_pos1 0 axis @vec height @height radius [jing_he*0.8]    endcommandend@ball_range

cmat default model rrlinear method deformability emod 100e6 kratio 1.5 property fric 0.5

ball generate radius 0.03 0.09 number 2000  tries 20000000 range geometry ballrange count 1ball attribute density 2.7e3 damp 0.7 set gravity 9.8 solvesave sample

生成后如图：

单向流固耦合——模拟弯曲河床砂石运移的图1

2 流场生成

这里使用我们之前的小程序 PFC单向流固耦合——模拟颗粒落入流动的水中生成流体网格。

参数为：

单向流固耦合——模拟弯曲河床砂石运移的图2

记得将其放置到项目文件下。

导入网格后我们需要对流场进行设置，这里我们假设流速是均匀的，实际上的河床截面上的流速是上大下小的应该。因为我们河床圆弧的中心在原点，所以我们直接根据流体网格的位置取切向作为流速方向即可。

restore sampleconfigure cfdball attribute displacement multiply 0cfd read nodes Node.datcfd read elements Elem.datcfd buoyancy on

[sudu=5]def addliusu    loop foreach local ele element.cfd.list        vec=math.unit(vector(element.cfd.pos.y(ele),(-1)*element.cfd.pos.x(ele),0))        element.cfd.vel(ele)=(-1)*vec*sudu    endloop end@addliusuelement cfd attribute density 1000.0element cfd attribute viscosity 1.5

set mech age 0[baocunpinlv=0.1][time_record=mech.age-1][count=0]def savefile        if mech.age-time_record >= baocunpinlv then        filename=string.build("jieguo%1",count)        command            save @filename        endcommand        time_record=mech.age        count +=1    endif    endset fish callback -1.0 @savefilesolve time 10