FLAC3D前后处理及三维模型构建方法

1、 岩土工程研究所Geotechnical Institute石石 崇崇河海大学岩土所河海大学岩土025-边坡岩石力学与工程边坡岩石力学与工程FLAC3D前后处理和三维模型构建方法Flac3d前后处理和三维模型构建方法主要认识和结论 FLAC3D计算结果的后处理 FLAC3D计算中的应力场FLAC3D三维模型的快速构建Sufer的三维出图和数据导出目录CAD的前处理和数据导出 在CAD中选择建模区域命令: Trim 对等高线截断 Move 平移,减小坐标起始值 导出为R12格式的dxf文件CAD的数据处理的数据处理 借助第三方软件GID导入CAD生成的dxf文件,将CA

2、D中的等高线信息读出. GID可以作为CAD和ANSYS的过渡平台.其可转换DXF-IGES格式. 在GID软件中显示等高线上的点信息,将等高线的点列为文本输出,作为sufer三维可视化出图的高程点数据. sufer出图只需要三维的高程点数据,格式写成如下: X1,Y1,Z1 X2,Y2,Z2 也可利用其它软件或编程读出CAD中的等高线数据信息 Flac3d前后处理和三维模型构建方法主要认识和结论 FLAC3D计算结果的后处理 FLAC3D计算中的应力场FLAC3D三维模型的快速构建Sufer的三维出图和数据导出目录CAD的前处理和数据导出 在滤波器中过滤CAD中导出的局部失真地形点. 设置z

3、方向上某一高程区域内的地形点为有效地形点. 如模型范围实际地形高程在5002000以内,则输入 z 2000 在插值方法中选择合适的插值方法.注意sufer只能提供栅格化插值方法.如要插值不规则网格,需自己编程实现. 设置插值区域:X方向 Min Max Y方向 Min Max Sufer可以生成多种效果图形.对于三维地形的显示,主要采用线框图和表面图.其中表面图的效果最好,其插值效果受插值点个数和插值方法制约. 地图-三维旋转可以旋转查看地形,设置光照,网格可以显示不同的效果.注意在sufer插值后,生成一个后缀名为grd的文件,此为存贮栅格化插值点三维信息的数据文件. 用用sufer打开打

4、开生成的生成的dgx.grd文件文件,选择菜单选择菜单:文件文件-文件另存为文件另存为在保存类型中的在保存类型中的下拉菜单中悬着下拉菜单中悬着ASCII XYZ *.dat格式格式,即为即为sufer生成的插值点数生成的插值点数据信息据信息.选择选择-选项选项-网格网格信息信息,可以查看插可以查看插值的精度和网格值的精度和网格间距间距.flac3d执行命令流command.txt，即可读入sufer生成的地表数据文件,投影生成flac3d，ansys和3dec的模型的命令流。 公路公路1875mnAutoCAD在DWG文件中以实体名存储各种点、线、面、实体。实体名是AutoCAD中专门设置的一

5、种特殊数据类型，是赋予图形实体的数字标号。通过AutoLisp语言开发的程序可以访问各种点、线、面、实体信息，例如边坡中的地形线，开挖面，地下厂房，引水洞、平硐等 n现阶段常用的插值方法有：Kriging法、Delaunay法、Polymesh法、Shepard法、Multiquadric法等等。滑动最小二乘法（Moving Least Squares Method，简称MLSM）是一种精度较高的方法，其克服了经典最小二乘法在拟合中的不足，实质为方差泛函极小化序列与最小二乘法。nSufer中都是栅格插值方法, 所以如果要按给定任意网格插值的话,需要自己编程.Flac3d前后处理和三维模型构建方

6、法主要认识和结论 FLAC3D计算结果的后处理 FLAC3D计算中的应力场FLAC3D三维模型的快速构建Sufer的三维出图和数据导出目录CAD的前处理和数据导出导入sufer插值的栅格点,在cad中按精度需要建立多条多义线polyline.类似与ansys中的askin命令的做法,在CAD也可以生成由多条多义线过渡圆滑的nurbs曲面. CAD操作命令Loft 按顺序选择几条多义线-回车生成面域. 在CAD中建立底面域,由封闭的多义线围成面,利用REGION生成面域,采用EXTRUDE命令按Z方向拉伸成体.EXTRUDE命令按Z方向拉伸成体.让其延伸至实际地面以上一定距离.为后面的CAD中布

7、尔切割操作做准备.利用slice命令中的曲面切割功能,将底面拉伸的块体沿实际地貌切割为两个块体. 同理断层和地层的切割也采用slice中的三点法或曲面法.采用slice切割出的地层与地貌交线利用CAD中的EXPORT命令将生成的实体SOLID模型输出为SAT格式的数据文件. ANSYS可以直接读入CAD导出的SAT模型. 在ANSYS进行网格剖分,建议全部剖分为四面体,可以设置线的剖分段数来控制局部和整体的网格质量. 由于FLAC3D中最终计算的单元均为四面体单元,因此采用四面体剖分精度影响不大.nFLAC3D数值模型分解后的最小单元类型为Brick（六面体）、Wedge（三棱柱体）、Pyra

8、mid（金字塔体）和Tetrahedral（四面体）单元。 在在ansys选择选择菜单菜单:FILE-Read Input from 选择选择APDL命令命令流文件流文件(在仿真论在仿真论坛可以下载此程坛可以下载此程序序)执行后在电脑的执行后在电脑的D:下生成下生成ansys单元节点信息的单元节点信息的两个文件两个文件01_node.dat02_ele.dat将将D:下生成下生成ansys单元节点单元节点信息的两个文件信息的两个文件01_node.dat02_ele.dat拷至程序目录拷至程序目录下下,运行运行flac_6weishu.exe后后,程序生程序生成了成了flac3d可可以直接读入

9、的以直接读入的模型数据模型数据.打开打开flac3d,读入生成的模读入生成的模型数据型数据. Flacmodel_haitang.Flac3D 读入后读入后,保存保存为为sav文件文件.结合CAD的extrude,slice,union,subtract,intersect等三维布尔操作功能,编制VBA,AUTOLISP程序可以快速构建更为复杂的实体模型.保存为sat格式后,可以直接被ansys读入.Flac3d前后处理和三维模型构建方法主要认识和结论 FLAC3D计算结果的后处理 FLAC3D计算中的应力场FLAC3D三维模型的快速构建Sufer的三维出图和数据导出目录CAD的前处理和数据导

10、出 岩体地应力是岩体工程最基本也是最重要岩体地应力是岩体工程最基本也是最重要的工程荷载之一。的工程荷载之一。 因此，初始应力分布应该作为设计和分析因此，初始应力分布应该作为设计和分析的基础，相应地，也应当作为数值计算的初始的基础，相应地，也应当作为数值计算的初始条件。条件。 初始应力场的重要性初始应力场的重要性ini sxx=-5e6 syy=-1e7 szz=-2e6Ini sxx s0 grad gx gy gz The components 11 (or xx), 22 (or yy) and 33 (or zz) are set to compressive stresses of -

11、5.0e6, 1.0e7 and -2 .0e6, respectively, throughout the grid. The INI command sets all stresses to the given values, respectively, throughout the grid. but there is no guarantee that the stresses will be in equilibrium.ngen zone brick size 6 6 6nmodel elasnini sxx=-5e6 syy=-1e7 szz=-2e7napply sxx=-5e

12、6 range x -0.1 0.1napply sxx=-5e6 range x 5.9 6.1napply syy=-1e7 range y -0.1 0.1napply syy=-1e7 range y 5.9 6.1napply szz=-2e7 range z -0.1 0.1napply szz=-2e7 range z 5.9 6.1Example 3.20 Initial and boundary stresses in equilibrium赋予初始应力ngen zone brick size 6 6 6nmodel elasnprop bulk 5e9 shear 3e9n

13、ini sxx=-5e6npausensolve只用ini，所有单元将有sxx应力，但solve变为0，因为没有约束，力引起位移最后平衡.所以必须加上应力边界条件或者位移边界条件加上应力边界条件ngen zone brick size 6 6 6nmodel elasnprop bulk 5e9 shear 3e9 nini sxx=-5e6napply sxx=-5e6 range x -0.1 0.1napply sxx=-5e6 range x 5.9 6.1npausenpl con sxxnsolvenpl con sxxn注意是约束x方向-变形参数并不影响应力场，它只影响位移场。如

14、果考虑重力的时候，变形参数影响泊松比，从而影响侧压力系数，才会影响水平的应力场分布。表明所有单元只加上SXX= -5e6Solve 是平衡的，即应力不发生改变了。改变bulk，shear和边界条件，solve后应力不发生变化。由于施加的初始应力场和应力边界条件产生了平衡，所以solve不产生cyc应力循环。只有不平衡力大于容许值时，才会进行cyc。设置不同的变形参数并不影响位移结果，因为系统平衡了。如果人为设置执行计算流程（如cyc 10000），由于没有位移边界条件约束，模型则会慢慢随时间产生位移，但应力场是不会发生变化的。因此，intial应力得到的位移解相比较于无intial应力的位移解

15、更符合实际地层条件（即认为是0）Cyc 20000或 step 20000ngen zone brick size 6 6 6nmodel elasnprop bulk 5e2 shear 3e2nini sxx=-5.0e6nfix x range x -0.1 0.1 nfix x range x 5.9 6.1 npausenpl con sxxnsolvenpl con sxxn注意是约束注意是约束x方向方向加上位移边界条件表明所有单元只加上SXX= -5e6Solve 是平衡的，即应力不发生改变了。改变bulk，shear和边界条件，solve后应力不发生变化。施加应力场后，如果没有

16、应力和位移约束，模型会沿着临空面产生变形。如果施加了位移边界条件约束了模型的变形，在边界上则产生了反向的平衡力，从而产生了平衡，所以solve不产生cyc应力循环。只有不平衡力大于容许值时，才会进行cyc。设置不同的变形参数并不影响位移结果，因为系统平衡了。如果人为设置执行计算流程（如cyc 10000），由于位移边界条件约束，边界上位移为0，内部则会慢慢随时间产生位移，但应力场是不会发生变化的。nnewn;Example 3.24 Initial stress state for an irregular free surfacengen zone brick size 15 15 10 p

17、0 0,0,0 edge=100.0nmodel elasticnprop shear 3e8 bulk 5e8ndef mountainngp = gp_headnloop while gp # nullnzz = sqrt(gp_xpos(gp)2 + gp_ypos(gp)2)ndz = 0.06 * sin(0.2 * zz + 100.0) ; Sum Fourier terms forndz = dz + 0.06 * sin(0.22 * zz - 20.3) ; quasi-random surfacendz = dz - 0.04 * sin(0.33 * zz + 33.3

18、) ; topology.ngp_zpos(gp) = 0.5 * gp_zpos(gp) * (1.0 + dz)ngp = gp_next(gp)nend_loopnendnmountainnfix x range x -.1 .1nfix x range x 99.9 100.1nfix y range y -.1 .1nfix y range y 99.9 100.1nfix z range z -.1 .1模型高度约为模型高度约为50m模型长宽约为模型长宽约为100m仅考虑自重，算至平衡仅考虑自重，算至平衡 (模型高度约为模型高度约为50m）泊松比泊松比=0.25 =0.33333

19、模型最大高程模型最大高程50m左右左右自重来看，基底自重来看，基底szz=2000*-10*50 = -1.0e6 sxx = syy = 0.333*-1.0e6=-0.333e6model elasticprop bulk 5e8 shear 3e8set grav 0,0,-10ini density=2000solve采用自重来计算初始应力场采用自重来计算初始应力场仅考虑自重，算至平衡仅考虑自重，算至平衡基底基底 szz=2000*-10*50 = -1.0e6 sxx = syy = 0.333*-1.0e6=-0.333e6prop bulk 5e8 shear 3e8prop b

20、ulk 5e11 shear 3e11此时此时 泊松比不变泊松比不变侧压力系数也不变侧压力系数也不变n按照自重求初始地应力，形成的应力场表现为竖直向=自重应力，水平向=自重应力侧压力系数k0。k0=v/(1-v). K0在这个方法中显然是和泊松比有关,然而土的泊松比变化不大,在0.25-0.35左右, 由此求得的k0在0.3-0.6之间,和phi在45-27度之间. 在上面,泊松比是通过语句prop bulk 5e9 shear 3e9间接来定义的.n这种方法会产生很大的位移,当然FLAC有一个投机取巧的方法来减小位移: 就是在求初始应力的时候,把材料设为弹性,并且用很大的弹性模量确(比如上面

21、的model elas; prop bulk 5e11 shear 3e11).修改变形参数不影响初始应力场的分布，只影响位移场的大小。n在此之后,初始化位移为0,并把材料的参数恢复正常. 但是这种方法存在两个问题：n第一是,由此求得初始侧压力系并不和实际一致,因此初始地应力也和实际不一致. 在这一点上用第二种方法可以得到比较正确的初始应力场. 或者也可以通过正确的k0反算一个不切实际的泊松比v代入来求.n第二也是由于设置材料为弹性的,所以是不会有破坏产生,不管初始产生的剪应力有多大,然而在材料设为弹塑性时,如果剪应力过大,材料就会在初始应力(没有其它外力作用)自动破坏,从而产生很大的变形,

22、因此有时候第一步计算时就会产生很大的变形,很有可能就是这个原因. 采用弹性材料计算的应力场在地形比较陡峭的地方有比较大的Mises应力，导致当采用弹塑性材料时破坏而不收敛（特别是有软弱夹层的地方）。那么，土体（土层）在原位应力下，土体基本上处于弹性状态。因此用弹性的方法估计初始应力是否可以满足计算的要求？n并不全然,在现实中土也有处于屈服状态的. 象边坡的潜在滑动面应力状态, 就很有可能处于殘余应力状态.为证明这一点,可以在FLAC建一个边坡模型,在线弹性模型下求得平衡, 然后再把弹塑性模型换回,在不加任何外载的情况下再求解一次,我相信会有很大的位移产生,并且会出性塑性屈服点.并可以预测,边坡

23、的安全系数越小,出现的塑性屈服点越多.那么，采用塑性的本构模型来算自重应力场，是否更为贴切呢？ini szz=-1.00e6 (grad 0,0,0.020e6)ini sxx=-0.36e6 (grad 0,0,0.072e6) ini syy=-0.36e6 (grad 0,0,0.072e6)solve手动构建初始应力场手动构建初始应力场利用利用ini设置相同于自重形成的初始应力场设置相同于自重形成的初始应力场自重应力场自重应力场 szz=1.0e6; sxx=syy=0.36e6 从公式上看从公式上看 地表地表szz = -1.0e6 + 50 * 0.02e6 = -0.0e6 sx

24、x = -0.36e6 + 50 * 0.072e6 = -0.0e6The INI command sets all stresses to the given values, respectively, throughout the grid. but there is no guarantee that the stresses will be in equilibrium.利用自重应力计算的初始应力场利用自重应力计算的初始应力场. szz=1.0e6; sxx=syy=0.36e6 利用利用ini设置的初始应力场设置的初始应力场. 未计算未计算solve。系统未平衡时。系统未平衡时in

25、i szz=-1.00e6 (grad 0,0,0.020e6) 整体规律与自重应力场相同整体规律与自重应力场相同ini sxx=-0.36e6 (grad 0,0,0.0072e6) ini syy=-0.36e6 (grad 0,0,0.0072e6)利用利用ini设置的初始应力场设置的初始应力场. 计算计算solve至系统平衡时，发现应力场降低了一个至系统平衡时，发现应力场降低了一个数量级，说明位移边界条件并不能使系统平衡，数量级，说明位移边界条件并不能使系统平衡，如何设置才能平衡呢？如何设置才能平衡呢？利用利用ini设置的初始应力场设置的初始应力场. 未计算未计算solve。系统未平衡

26、时。系统未平衡时ini szz=-1.00e6 (grad 0,0,0.020e6) 整体规律与自重应力场相同整体规律与自重应力场相同ini sxx=-0.36e6 (grad 0,0,0.0072e6) ini syy=-0.36e6 (grad 0,0,0.0072e6)以上分析可知，单纯用自重来计算应力场，耗时较长。如果仅通过以上分析可知，单纯用自重来计算应力场，耗时较长。如果仅通过ini手动设置应力场，手动设置应力场，有时并不能使系统平衡，这说明加入人为得到的初始地应有时并不能使系统平衡，这说明加入人为得到的初始地应initial，在弹性情况下是不，在弹性情况下是不会对初始应力造成影响

27、的。那会对初始应力造成影响的。那如何设置才能平衡呢？如何设置才能平衡呢？可以想像，可以想像，ini是人为设置的一个单元应力，其可以加快平衡迭代的过程。是人为设置的一个单元应力，其可以加快平衡迭代的过程。在计算初始在计算初始应力场的时候，如果将岩体设定为弹性的话，在给定的重力与约束条件下，初始应力应力场的时候，如果将岩体设定为弹性的话，在给定的重力与约束条件下，初始应力是不会对最终的应力结果造成影响的。是不会对最终的应力结果造成影响的。ini是模型的内力，是模型的内力，set grav是外力施加对模型是外力施加对模型的作用力，没有外力怎么有内力？这样做其实只是减少程序的计算步骤的作用力，没有外力

28、怎么有内力？这样做其实只是减少程序的计算步骤. 这样设置以这样设置以后后 整个模型就已经是平衡了。同时设置重力的原因不是为了平衡减少程序计算步骤，整个模型就已经是平衡了。同时设置重力的原因不是为了平衡减少程序计算步骤，而是因为如果未初始化应力的话，在而是因为如果未初始化应力的话，在set grav的时候，模型会在自重应力下发生应力的时候，模型会在自重应力下发生应力变形。所以，初始化应力的原因是防止模型在自重应力下发生变形变形。所以，初始化应力的原因是防止模型在自重应力下发生变形 一、施加自重，施加位移边界条件一、施加自重，施加位移边界条件model elasticprop bulk 5e11

29、shear 3e11set grav 0,0,-10ini density=2000ini szz=-1.00e6 (grad 0,0,0.020e6)ini sxx=-0.36e6 (grad 0,0,0.0072e6) ini syy=-0.36e6 (grad 0,0,0.0072e6)fix x range x -.1 .1fix x range x 99.9 100.1fix y range y -.1 .1fix y range y 99.9 100.1fix x y z range z -.1 .1Solve 需需527步平衡步平衡一、施加自重，施加应力边界条件一、施加自重，施加

30、应力边界条件model elasticprop bulk 5e11 shear 3e11set grav 0,0,-10ini density=2000ini szz=-1.00e6 (grad 0,0,0.020e6)ini sxx=-0.36e6 (grad 0,0,0.0072e6) ini syy=-0.36e6 (grad 0,0,0.0072e6)apply sxx -0.36e6 (grad 0,0,0.0072e6)apply syy -0.36e6 (grad 0,0,0.0072e6)fix x y z range z -.1 .1Solve 需需875步平衡步平衡一、施加

31、自重，施加应力边界条件一、施加自重，施加应力边界条件model elasticprop bulk 5e11 shear 3e11set grav 0,0,-10ini density=2000ini szz=-1.00e6 (grad 0,0,0.020e6)ini sxx=-0.36e6 (grad 0,0,0.0072e6) ini syy=-0.36e6 (grad 0,0,0.0072e6)apply sxx -0.36e6 (grad 0,0,0.0072e6)apply syy -0.36e6 (grad 0,0,0.0072e6)fix x y z range z -.1 .1S

32、olve 需需875步平衡步平衡一、施加自重，施加位移边界条件一、施加自重，施加位移边界条件model elasticprop bulk 5e11 shear 3e11set grav 0,0,-10ini density=2000ini szz=-1.00e6 (grad 0,0,0.020e6)ini sxx=-0.36e6 (grad 0,0,0.0072e6) ini syy=-0.36e6 (grad 0,0,0.0072e6)fix x range x -.1 .1fix x range x 99.9 100.1fix y range y -.1 .1fix y range y 9

33、9.9 100.1fix x y z range z -.1 .1Solve 需需527步平衡步平衡在某些构造区域，很可能水平应力远大于竖直应力，此时不能用自重来计算在某些构造区域，很可能水平应力远大于竖直应力，此时不能用自重来计算初始应力场。此时采用初始应力场。此时采用ini来影响水平应力场的分布。来影响水平应力场的分布。There is no simple way to deduce an equilibrium stress distribution for this grid: it must be modeled as a boundary-value problem. Howeve

34、r, we may insert initial stresses in order to speed up the convergence and influence the final stress distribution. For example, if we know that there is a high horizontal in-situ stress, with only a small decrease in stress near the surface, we can initialize xx and yy to -0.2e6 at the approximate

35、location of the surface, increasing to -4.0e6 at the bottom. The vertical stress, zz, can be set to correspond to the average overburden.set grav 0,0,-10ini density=2000ini szz=-1.00e6 (grad 0,0,0.020e6)ini sxx=-4.00e6 (grad 0,0,0.076e6) ini syy=-4.00e6 (grad 0,0,0.076e6)fix x range x -.1 .1fix x ra

36、nge x 99.9 100.1fix y range y -.1 .1fix y range y 99.9 100.1fix x y z range z -.1 .1Solve 需需402步平衡步平衡ini sxx=-4.00e6 (grad 0,0,0.076e6) ini syy=-4.00e6 (grad 0,0,0.076e6)地表地表sxx = -4.0e6 + 50 * = -0.2e6所以所以 = 0.076e6由此可见，由此可见，ini是人为设置的一个单元应力，不仅可以加快平衡迭代的过程，而且可以采用是人为设置的一个单元应力，不仅可以加快平衡迭代的过程，而且可以采用ini来影

37、响水平应力场的分布。同时设置的来影响水平应力场的分布。同时设置的ini水平初始应力是不会对竖直应力场造成影响的。一水平初始应力是不会对竖直应力场造成影响的。一般而言这样是为了计算比较快般而言这样是为了计算比较快,因为这样设置以后因为这样设置以后,此时几乎已经平衡此时几乎已经平衡,了了,计算步较少了计算步较少了,变形变形也就很少了嘛也就很少了嘛,当然最好再在以后计算时初始化位移为当然最好再在以后计算时初始化位移为0.当然还有另一种方法就是开始时设成当然还有另一种方法就是开始时设成弹性的弹性的(model elas)（或者还是（或者还是mohr的但参数给的要大一些，目的是防止破坏），直接运行的但参

38、数给的要大一些，目的是防止破坏），直接运行到平衡达来获取初始应力值，在真正计算时再改成到平衡达来获取初始应力值，在真正计算时再改成mohr并将参数改为真实的数值，这种方并将参数改为真实的数值，这种方法更适合于单元不太多，因为比较简单但如果单元数目很多的话，那么此种方法速度会很法更适合于单元不太多，因为比较简单但如果单元数目很多的话，那么此种方法速度会很慢慢的的n1. Do not initialize stresses; allow gravity to compact the layer.n2. Initialize horizontal stress only, not vertical

39、stress.n3. Impose constant stress at the lateral boundaries rather than zero horizontal displacement.n4. Remove irregular overburden from initial grid of uniform thickness.n5. Allow plastic flow to occur, thus removing stress concentrations.n6. Build up the profile layer by layer; equilibrate each l

40、ayer.n7.初始地应力的问题，应该说是一个困难的问题，并没有一个完整的解决方案。在有限元工具中，比如ansys或者abaqus，都可以先固定边界，然后施加重力得到一个弹性应力场，然后再将这个应力场作为初始应力场导入和重力平衡，的确可以不产生位移。但关键的问题是，采用弹性材料计算的应力场在地形比较陡峭的地方有比较大的米赛斯应力，导致当采用弹塑性材料时破坏而不收敛（特别是有软弱夹层的地方）。往往这些区域都不大，在靠近地表比较陡峭的地方。虽然存在这些缺点，但仍不失为一个好方法。采用ini条件加线性的初始应力，这个方法flac有，abaqus也有，但这个方法的缺点不言而喻，不合适地形起伏的情况，但

41、最起码他可以控制Ko，这点又是前面这方法所欠缺的。所以，这两招都有其合适的场合，关键是看工程的实际作用范围了，在工程的影响范围之外的地应力影响比较小，是否可以比较粗糙一些。n1. 利用重力来计算竖直向应力场利用重力来计算竖直向应力场n2. 人为施加水平向应力场人为施加水平向应力场,加快收敛速度加快收敛速度,同时可以得到实际的同时可以得到实际的水平应力场分布水平应力场分布n3. 采用塑性本构计算初始应力场采用塑性本构计算初始应力场,允许塑性流动发生允许塑性流动发生,这样可以这样可以减小局部应力集中的现象减小局部应力集中的现象 .n4. 等比例提高等比例提高bulk和和shear参数参数,是不会改

42、变泊松比的是不会改变泊松比的.这样在这样在不改变应力场的前提下不改变应力场的前提下,也减小了自重产生的位移场也减小了自重产生的位移场.(这也是自这也是自重产生的应力场在后面的计算中是要清零的原因重产生的应力场在后面的计算中是要清零的原因.).n5.如有原位测点的地应力资料如有原位测点的地应力资料.分段计算出不同区位的梯度施分段计算出不同区位的梯度施加力和梯度值加力和梯度值. 由于由于flac以原点来计算梯度力以原点来计算梯度力 ini sxx s0 (grad 0,0,gz) rang z z1 z2测点应力场测点应力场 h1处为处为s1; h2处为处为s2 从公式上看从公式上看 sxx/h1

43、 = s1 = s0 + 0*x1 + 0*y1 + gz*z1 sxx/h2 = s2 = s0 + 0*x2 + 0*y2 + gz*z2 河谷边坡中的地应力场是一个形如岩体开挖河谷边坡中的地应力场是一个形如岩体开挖以后的二次应力场，不适宜用构造应力分量与以后的二次应力场，不适宜用构造应力分量与自重应力分量的叠加；自重应力分量的叠加； 理论上，河谷发育是在构造运动格局基本确理论上，河谷发育是在构造运动格局基本确定以后的地表地质作用；模拟河谷侵蚀下切来定以后的地表地质作用；模拟河谷侵蚀下切来获得局部地应力场与真实地应力场更接近。获得局部地应力场与真实地应力场更接近。河谷地应力场模拟思想河谷地

44、应力场模拟思想313222113213河谷地应力场的基本特征河谷地应力场的基本特征n河谷地构造应力是一个非常不均匀的、受河谷发育改造的局部地应力场河谷地应力场模拟的基本考虑河谷地应力场模拟的基本考虑x坝址区两岸为单斜构造，右岸地形呈现5级阶面，因此，将河谷侵蚀下切假设为5期1300m高程 坝址区的构造应力主方向假设五期河谷侵蚀下切假设五期河谷侵蚀下切 地应力场模拟结果 正确性验证谷底高应力区（岩芯饼化）应力松弛带Flac3d前后处理和三维模型构建方法主要认识和结论 FLAC3D计算结果的后处理 FLAC3D计算中的应力场FLAC3D三维模型的快速构建Sufer的三维出图和数据导出目录CAD的前

45、处理和数据导出nFlac3d只能出云图只能出云图,不能出等值线图不能出等值线图.而工程中大多而工程中大多要求出等值线图要求出等值线图.nFlac3d三维截剖面出图效果不好三维截剖面出图效果不好,不能截几个剖面一不能截几个剖面一起出图起出图n所以所以,有必要导出其计算的应力位移信息有必要导出其计算的应力位移信息,利用利用sufer或或tecplot出图出图.nSufer:n优点优点:要求数据结构简单要求数据结构简单,可以采用白化功能进行局部区域的等值线出图可以采用白化功能进行局部区域的等值线出图.n缺点缺点:需要给定边界需要给定边界(bln文件文件),不能出三维等值线不能出三维等值线.TECPL

46、OT:n优点优点:不需边界不需边界bln文件文件.可以出三维等值线图可以出三维等值线图.截剖面方便截剖面方便,可同时多剖面出图可同时多剖面出图n缺点缺点:数据结构复杂数据结构复杂,不能隐藏部分单元出图不能隐藏部分单元出图. DEFINE GetGPDispX P1_p=gp_head loop n (1,P_total) node_id=gp_id(p1_p) message(1)= string(gp_xpos(P1_p) message(1)=message(1)+,+string(gp_ypos(P1_p) message(1)=message(1)+,+string(gp_xdisp(

47、P1_p)*1000) message(1)=message(1)+,+string(gp_ydisp(P1_p)*1000) message(1)=message(1)+,+string(gp_zdisp(P1_p)*1000) status = write(message,1) P1_p=gp_next(P1_p) endloopENDGetGPDispX采用采用flacflac内嵌的内嵌的fishfish语言进行所有节点的位移值提取语言进行所有节点的位移值提取def initializationIO_READ = 0IO_WRITE = 1IO_FISH = 0IO_ASCII = 1S

48、ufer_dispX =string( Sufer_dispX.txt )array message(1)endinitializationstatus = open(Sufer_dispX,IO_WRITE,IO_ASCII) if status = 0 then GetGPDispX status = close ii = out(Successfully Write Data Into File + Sufer_dispX)else ii=out(Open File Error! Status = + string(status)endif采用采用flacflac内嵌的内嵌的fishfish语言进行所有节点的位移值提取语言进行所有节点的位移值提取也可采用也可采用set log onset log on的方式记录位移的方式记录位移. .但是但是loglog文件中文件中也包含了其它信息也包含了其它信息. .建议采用文本通道加格式输出的建议采用文本通道加格式输出的方式方式, ,比较方便比较方便. .; originl by dynamaxsimwe, modified by HaitangHohai; 18/12/