FLAC3D网上学习资料整理3

FLAC3D网上学习资料整理3

网上FLAC3D学习笔记

15.建立interface的几种方法 首先给出一个简单的例子 new

gen zone brick p0 0 0 0 p1 111.4 0 0 p2 0 10 0 p3 0 0 45 size 11 1 5 group 1

gen zone brick p0 0 0 45 p1 111.4 0 45 p2 0 10 45 p3 0 0 90 p4 111.4 10 45 &p5 0 10 90 p6 47 0 90 p7 47 10 90 size 11 1 9 group 3

gen zone wedge p0 66.4 0 90.1 p1 47 0 90.1 p2 66.4 10 90.1 p3 111.4 0 45.1 &p4 47 10 90.1 p5 111.4 10 45.1 size 5 1 9 group 4

;interface 1 face rang plane norm 45,0,64.4 origin 47 5 90 interface 1 face range plane norm 45,0,64.4 origin 47 5 90 dist 0.1

ini z add -0.1 range group 4

gen zone brick p0 111.4 0 0 p1 128 0 0 p2 111.4 10 0 p3 111.4 0 45 &

size 2 1 5 group 2 model mohr

pro bulk=1e8 she=3e7 dens 2e3 coh 12e3 & friction 20 dilation 0. tension 1.0e10 range group 1 pro bulk=1e8 she=3e7 dens 2e3 coh 12e3 & friction 20 dilation 0. tension 1.0e10 range group 2 int 1 prop ks 2e12 kn 2e12 fric 36

pro bulk=1e8 she=3e7 dens 2e3 coh 12e3 & friction 20 dilation 0. tension 1.0e10 range group 3 pro bulk=1e8 she=3e7 dens 2e3 coh 12e3 & friction 20 dilation 0. tension 1.0e10 range group 4 fix z range z -0.1 0.1 fix x range x -0.1 0.1 fix x range x 127.9 128.1 fix y

set gravity 0 0 -10 step 1000 ini state 0

ini xdis 0 ydis 0 zdis 0 ;ini xvel 0 yvel 0 zvel 0 ;set large step 1000

;solve fos file cable.sav associated 修改上面蓝色字体部分

第一种建立interface的方法（书上的移来移去式）：

gen zone brick p0 0 0 0 p1 111.4 0 0 p2 0 10 0 p3 0 0 45 size 11 1 5 group 1

1

pl int yel sket

通过计算得到的xdisp云图

gen zone brick p0 0 0 45 p1 111.4 0 45 p2 0 10 45 p3 0 0 90 p4 111.4 10 45 &p5 0 10 90 p6 47 0 90 p7 47 10 90 size 11 1 9 group 3

gen zone wedge p0 66.4 0 90.1 p1 47 0 90.1 p2 66.4 10 90.1 p3 111.4 0 45.1 &p4 47 10 90.1 p5 111.4 10 45.1 size 5 1 9 group 4

;interface 1 face rang plane norm 45,0,64.4 origin 47 5 90 interface 1 face range plane norm 45,0,64.4 origin 47 5 90 dist=0.1

ini z add -0.1 range group 4

gen zone brick p0 111.4 0 0 p1 128 0 0 p2 111.4 10 0 p3 111.4 0 45 &

size 2 1 5 group 2

FLAC3D网上学习资料整理3

问：为什么第一种方法要加 dist 0.1？ 加入节理以后，这个例子就不能计算fos 了？

回答：你的第2种，把group2和group4删掉，这样造成的新的模型group2和group4有重复节点，这样是不正确的。方法是在一个实体的外表面全部都设置了接触面，就是可以将两个体完全分开。而你的这个接触问题，在坡脚的位置很难处理。只用把斜坡上的部分删掉，然后导入就可以了。

人为的把group1和group2之间的交界面上多出了节点，下图是两个group的节点图，是分开显示的，可以看到两个共同面上界面标号不同，这样是不对的。计算结果也同样是错误的。

behave elastically during the development of the initial in-situ stress state.* This prevents any plastic yield during the initial loading phase of the analysis. Gen zone brick size 6 8 8 Mode mohr

Prop bulk 1e8 shear 0.3e8 fric 35

Prop cohesion 1e10 tens 1e10 ;注意在此这个值给的很大。 Init dens 1000 Set gravity 0 0 -10 Fix x range x -0.1 0.1 Fix x range x 5.9 6.1 Fix y range y -0.1 0.1 Fix y range y 7.9 8.1 Fix z range z -0.1 0.1 set mech force=50 solve

;---------------------- mode null ---------------------

Prop coh 1e3 tens 1e3 ;改为正常值(在此例中我们故意给小值) Mode null range x 2 ,4 y 2 , 6 z 5, 10 Set large

16.call fishcall.fis总是会出错呢？老是给出error opening file，该怎么操作？

错误是：文件名为 fishcall.fis.txt(因为在记事本里这个txt的后缀是默认的）。修改方法是：将 fishcall.fis.txt里面的内容全部复制到新的记事本里，再按保存，保存时候将保存类型改为所有文件(*,*),这时你便可在文件名栏内输入 fishcall.fis，类型也就正确了 17.问题：

fix z range z -0.1 0.1

fix x range x -0.1 0.1 any x 1.9 2.1 any fix y range y -0.1 0.1 any y 1.9 2.1 any 和

fix z range z -0.1 0.1 fix x range x -0.1 0.1 fix x range 1.9 2.1 fix y range y -0.1 0.1 fix y range y 1.9 2.1

问：两者有区别吗？我感觉是一样的

回答：建一个模型，然后实验一下，然后用plot gpfixt看一下。 18. 得到初始应力的方法：

方法1、可以先给一些材料参数很大的值，进行初始求解，在计算之前再将材料参数设为正常值，即可。如在手册中给的第一个示例中就是这样做的。下面是例子，These are only initial values that are used during the development of gravitational stresses within the body. In effect, we are forcing the body to

2

Ini xdis 0 ydis 0 zdis 0 ;清零，不影响结果，为画图方便。 Setp 2000 return

方法2：得到初始应力的一个不错办法是设为弹性模型，然后赋予地勘信息，so以后就ok了，得到初始应力了，再把相应的位移清零，下面就可以正式进入所分析的内容了。如上例还可以先设为弹性 (MODEL elastic)，之后在开挖之前再设为 Mohr-Coulomb model， 就行。 19．关于初始应力的问题

Example 3.21 Initial stress state with gravitational gradient for example, a 20 m 20 m 20 m box of

homogeneous material at a depth of 200 m underground, with fixed base and stress boundaries onthe other sides. gen zone brick size 10,10,10 p1 20,0,0 p2 0,20,0 p3 0,0,20 model mohr

prop bulk 5e9 shear 3e9 fric 35 ini density 2500 set gravity 0,0,-10 fix x y z range z -0.1 0.1 ini szz = -5.0e6 grad 0,0,2.5e4 ini sxx = -2.5e6 grad 0,0,1.25e4 ini syy = -2.5e6 grad 0,0,1.25e4 apply szz = -4.5e6 range z 19.9 20.1 apply szz = -5.0e6 range z -0.1 0.1

apply sxx = -2.5e6 grad 0,0,1.25e4 range x -0.1 0.1 apply sxx = -2.5e6 grad 0,0,1.25e4 range x 19.9 20.1

FLAC3D网上学习资料整理3

apply syy = -2.5e6 grad 0,0,1.25e4 range y -0.1 0.1 apply syy = -2.5e6 grad 0,0,1.25e4 range y 19.9 20.1 它的意思是通过ini szz来赋予初始地应力,但是为何又有set gravity 0,0,-10?这样不是重复了吗?

ini是模型的内力，set grav是外力施加对模型的作用力，没有外力怎么有内力？一般而言，这样是为了计算比较快,因为这样设置以后,此时几乎已经平衡,了,计算步较少了,变形也就很少了嘛,当然最好再在以后计算时初始化位移为0,但初始化位移与否(在这种情况下),并不会明显影响计算结果的.当然还有另一种方法就是开始时设成弹性的(model elas)（或者还是mohr的但参数给的要大一些，目的是防止破坏），直接运行到平衡达来获取初始应力值，在真正计算时再改成mohr并将参数改为真实的数值，这种方法更适合于单元不太多，因为比较简单．但如果单元数目很多的话，那么此种方法速度会很慢的． 关于初始地应力及其算法，有两个疑问： 1）为何一定要让模型达到初始平衡再进行计算？

2.）那种改大参数，让模型在自重下平衡得到初始应力场的做法中，要将哪些参数改大呢？bul she coh ？还是改其他哪几个？ 回答：

1）这是模拟多数现场施工前的实际情况,并非总是如此,例如土钉墙施工,如果开挖后立即加土钉,就无须先平衡再加土钉再SOLVE.

2) 取决于具体模型,核心思想是用该模型在此阶段模拟弹性模型.对于MOHR-COULOMB模型, 这样做(将

COHESION,TENSION设置很大的值)与SOLVE ELASTIC是等价的.但对于别的模型,SOLVE ELASTIC并不适用,最好用改变参数的方法来加快计算速度和防止屈服.另外一点: 在FLAC/FLAC3D中, SET命令作用于整个网格(赋值全局变量),INI则初始化指定范围内网格点/域的偏移量.初始应力不仅仅是重力引起的。它还包括土压力，水浮力，上方建筑物的力，这些信息的勘测一般都不容易获得。通常可以把模型材料设为弹性，再把两个模量设到很大solve后再把模量和材料还原可得。具体可以看以前关于初始应力的帖子。不只是这样的，set grav不是直接将力仅仅加在底层上，而是作用在上面，如果没有初始应力的话，地层就会有相应的变形直到地层变形和产生的力与之平衡。所以应当通过initial来初始化应力，而这个力与density相匹配。我是做地铁开挖的如果没有初始化的话地层在开挖前就有了沉降。这与实际不符合。设置重力的原因不是为了平衡减少程序计算步骤，而是因为如果未初始化应力的话，在set grav的时候，模型会在自重应力下发生应力变形。所以，初始化应力的原因是防止模型在自重应力下发生变形（我们认为地层自重变形已经发生，当然欠固结等现象是另外一回事）。其实就是让模型发生变形，然后把位移清零也可以的。有种方法，不少人也在用，就是为了得到初始应力值，先将材料参数放大到很大，获得初始值以后再恢复到正常参数。

3

20. Tension-p shear-p tension-n shear-n表示的含义 在计算的循环里面,每个循环中,每个zone都依据failure criterion处于不同的状态,shear 和tension分别表示因受剪和受拉而处于塑性状态.n表示now,p表示previous,即分别表示在现在和以前的循环当中处于塑性状态. 21. 有一个例子：

gen zone cyl p0 0 0 0 p1 1 0 0 p2 0 2 0 p3 0 0 1 size 4 5 4 gen zone reflect norm 1,0,0 gen zone reflect norm 0,0,1 model mohr

prop bulk 1.19e10 shear 1.1e10 prop coh 2.72e5 fric 44 ten 2e5 fix x y z range y -.1 .1 fix x y z range y 1.9 2.1 ini yvel 1e-7 range y -.1 .1 ini yvel -1e-7 range y 1.9 2.1 hist gp ydisp 0,0,0 hist zone syy 0,1,0 hist zone syy 1,1,0

hist write 1 file 原点位移.dat hist write 2 file 中部单元应力.dat hist write 3 file 边缘单元应力.dat hist write 1 table 1 hist write 2 table 2 hist write 3 table 3

hist write 2 3 vs 1 file 单压.dat step 3000

结果在“原点位移.dat”文件中只有：“Step 1 Y-Displace” 而没有保存下数值，怎么回事？我怎么才能获得table中的数据？ 修改：

gen zone cyl p0 0 0 0 p1 1 0 0 p2 0 2 0 p3 0 0 1 size 4 5 4 gen zone reflect norm 1,0,0 gen zone reflect norm 0,0,1 model mohr

prop bulk 1.19e10 shear 1.1e10 prop coh 2.72e5 fric 44 ten 2e5 fix x y z range y -.1 .1 fix x y z range y 1.9 2.1 ini yvel 1e-7 range y -.1 .1 ini yvel -1e-7 range y 1.9 2.1 hist gp ydisp 0,0,0 hist zone syy 0,1,0 hist zone syy 1,1,0 step 3000

hist write 1 file 原点位移.dat

FLAC3D网上学习资料整理3

hist write 2 file 中部单元应力.dat hist write 3 file 边缘单元应力.dat hist write 1 table 1 hist write 2 table 2 hist write 3 table 3

hist write 2 3 vs 1 file 单压.dat

save a 'print tabel n' will solve your second qusetion 22.这个plot是什么意思？ hist n 1

hist gp ydisp 0,0,0 hist ax_str hist gp xdisp 1,1,0 step 1500

plot hist 2 vs -1 ;axial stress vs axial disp. plot hist 3 vs -1 ; circumferential disp. vs axial disp

plot 显示 hist 历史曲线 曲线2为横坐标，曲线一为纵坐标 hist n 1

hist gp ydisp 0,0,0 hist ax_str hist gp xdisp 1,1,0 step 1500

plot hist 2 vs -1 ;axial stress vs axial disp.

plot hist 3 vs -1 ; circumferential disp. vs axial disp

那2是不是代表hist gp ydisp 0,0,0，而 -1 指得是hist gp ydisp 0,0,0为什么不是1而是负值？ hist gp ydisp 0,0,0 hist ax_str hist gp xdisp 1,1,0 step 1500

plot hist 2 vs -1 ;axial stress vs axial disp. plot hist 3 vs -1 ; circumferential disp. vs axial disp 1) 在plot hist m vs n的形式里,m代表y轴,n代表x轴(不管m,n的正负);

2) "-"表示对其值作"mirror",如果大家使用AUTOCAD的话就明白"mirror"的含义了,比如上面的 hist 1 gp ydisp的值是0~100,那么vs -1就变成-100~0.以此类推.

3) 为什么要这么做呢?主要是为了符合岩土工程的习惯需要,例如桩载荷试验曲线,沉降曲线等都不是画在第一象限内. 23.关于water.fis的命令注释（设置地下水位）

4

1）water.fis的功能是设置地下水位面，在二维计算中，通过water table命令即可设置水位线，在三维计算中由于地下水位是面，因此需要作一些插值运算，以便形成水位面。这个FISH程序包含三个函数parm,aux1和aux2。parm用来设定水位的边值，即设置水位的区域，aux1插值计算等值线的点，aux2调用aux1的计算结果，使用FLAC3D命令形成水位面。 2）首先设置参数值，然后定义一个table. def parm nptab1=7 nprof=10 end parm

table 1 (12 0) (18.25 12.5) (25 17) (35 21) (50 24) (75 25) (150 25)

产生table1，包含7组数据点。

3) 在解释的两个FISH函数之前，对table,xtable,ytable作一些总的说明。water.fis程序反复使用了这三个变量，理解了这三个变量的含义，water.fis这段程序就解决了大半问题了。 Itasca的所有软件都使用了相同的数据结构储存数据，称之为链表结构，类似于关系型数据库中的数据链接方法。如果有数据结构的知识，或掌握关系型数据库的SQL语言，这三个变量的含义就不难理解了，只不过这个操作与真正的数据库操作相比，只有两个数据变量x,y。我用数据库的术语作解释：table就是产生一个“表”，xtable,ytable的作用就是产生一个新表(Create)以及对表中的数据进行取出(Select)、插入(Insert)、替换(Update)操作.

4）FISH函数aux1的作用是在table 1的基础上产生出另外个table，table的ID从11到17，每个table有10个数据集。alfa的计算方法是为适合(fit)曲面而采用的坐标点算法。 def aux1 loop n(1,nptab1) rr=xtable(1,n) zz=ytable(1,n) loop k(1,nprof) case_of k

alfa=0.5*(pi+.1)*float(k-2)/float(nprof-2) xtable(n+10,k)=rr*cos(alfa) ytable(n+10,k)=rr*sin(alfa) case 1

xtable(n+10,k)=rr ytable(n+10,k)=-40.0 case 2

xtable(n+10,k)=rr ytable(n+10,k)=0.0 end_case end_loop

FLAC3D网上学习资料整理3

end_loop end aux1

当执行这个FISH函数后，生成(nptab1+10,nprof)组数据，这些数据为Aux2调用。注意：本段程序产生的数据是为那个特定的曲面边坡制作的，不是一个general purpose的LIB 5）FISH函数aux2的作用是使用WATER TABLE face命令产生水

面，把table 1的ytable值作为z坐标，aux1产生的数据作为x,y坐标。 def aux2 command set grav 0 0 -10 water dens 1000 end_command loop n(1,nptab1-1) loop k (1,nprof-1) xx1=xtable(n+10,k) yy1=ytable(n+10,k) zz1=ytable(1,n) xx2=xtable(n+10,k+1) yy2=ytable(n+10,k+1) zz2=zz1

xx3=xtable(n+11,k+1) yy3=ytable(n+11,k+1) zz3=ytable(1,n+1) xx4=xtable(n+11,k) yy4=ytable(n+11,k) zz4=zz3 command

WATER WATER

TABLE TABLE

face face

xx1,yy1,zz1 xx1,yy1,zz1

xx2,yy2,zz2 xx3,yy3,zz3

xx3,yy3,zz3 xx4,yy4,zz4 end_command end_loop end_loop end aux2

24.产生结构元的几何形状(Geometry Creation) FLAC3D有6种类型的结构单元，它们分别是梁单元beamSELs,锚杆单元cableSELs,桩单元pileSELs,壳单元shellSELs,地质网格单元geogridSELs和衬砌单元linerSELs.每一种类型的结构由相应分量对象的集合来表示。例如锚杆结构(cable)由一系列的cableSELs组成，而衬砌结构(liner)则由

5

图中红色的数字1-3-4-2代表什么？是node？还是link？ sel pile id=1 begin=(10.0, 1.0, 0.0) end=(10.0, 1.0, -3.0) nseg=3

sel pile id=2 begin=(10.0, 3.0, 0.0) end=(10.0, 3.0, -3.0) nseg=3

sel pile id=3 begin=(10.0, 5.0, 0.0) end=(10.0, 5.0, -3.0) nseg=3

ID是为结构进行编号，CID是为结构的分量进行编号。比如下面的例子：

sel pile id=1 begin=(10.0, 1.0, 0.0) end=(10.0, 1.0, -10.0) nseg=5

sel pile id=2 begin=(10.0, 3.0, 0.0) end=(10.0, 3.0, -10.0) nseg=5

在这个例子中，设置了2根桩，第一根桩以id=1来表示，第二根桩以id=2来表示；第一根桩分为5段，则cid=1~5，第二根桩也分为5段，则cid=6~10。 可以用plot sel geom 命令来查看设置是否正确。 问题: 我设了三根桩，

CID:用来表示一个结构对象内部的分量，比如一根锚杆划分为3段，则这根锚杆的CID分别是1,2,3,用来表示这根锚杆的不同部分。

用Seltype命令和相应的关键字来定义结构单元:Seltype beam/cable/geogrid/liner/pile/shell,每一种支护类型的性质用SEL 命令的property关键字定义. History SEL keywords

ID：用来标识一个结构对象，比如第一根锚杆的ID＝1，第二根锚杆的ID=2，它是对整个结构进行编号；

一系列的linerSELs组成。每种结构元都由两个不同的辨识号码来表示，以便与其它结构元相区别。

生成node的时候，如果周围存在zone，会自动生成link，而且link和node位置总是重合的

本文来源：https://www.bwwdw.com/article/2fre.html