FLAC-讲义(含地下水的边坡例子).doc

FLAC 讲义一、什么是FLAC 1.1 FLAC之字义 F(Fast)L(Lagrangian)A(Analysis of)C(Continua). Lagrangian相对于Eulerian为每一时阶(timestep)之位 移在Lagrangian之公式中，需对网格之座标予以更新， 而Eulerian之公式则不予更新。1. 2 FLAC之运算流程 1.3 FLAC 基本单元1.4 分析模式大小与RAM之关系1.5 单位1.6 正负号方向 (1) 应力正号代表张力，负号代表压力 (2) 剪应力详见下图，图中所示剪应力为正号(3) 应变正的应变表示伸长，负的应变代表压缩(4) 剪应变剪应变的正负号与剪应力相同(5) 孔隙压力孔隙压力永远为正(6) 重力正号的重力物质往下拉，负号的重力将物质 往上提。二、FLAC内建之组合律FLAC内建之组合律有：1 空洞模式(null model) 使用于土壤被移除或开挖 2 弹性模式 3 塑性模式，包括a. Drucker Prager modelb. MohrCoulomb model c. ubiquitousjoint modeld. strainhardening/softening model e. bilinear strainhardening/softening modelf. doubleyield modelg modified camclay model此外，另有选购(option)模式，包括：1. 动力模式(Dynamic Option)2. 热力模式(Thermal Option)3. 潜变模式 (Creep Option) 使用者另可使用FISH语言去建构独特的组合律以符 合所需。三、FLAC以命令为输入语法请查阅相关手册四、FLAC程式之使用步骤 4.1 FLAC程式使用前准备步骤 步骤1：依比例画出所欲分析之资料于纸上画出地点之位置、地层资料、并简标示距离及深度资料。 步骤2：换算输入资料成同一单位将现有地层资料，如 Density,Bulk modulus,Youngs modulus, tension,cohesion, friction Angle 等资料，换算成同一单位。附注 ：需谨慎检查输入资料之单位，如因单位 不同而造成过大或过小的值，将会造成 FLAC无法计算，而产生ERROR讯息。 步骤3：应用公式简略计算应用公式或依据经验，简略算出FLAC 输 出资料的范围，以做为 Debug 及输出资 料分析时验证。步骤4：建立x,y座标与node i,j之关系 于图上距离及深度之关系，建立x,y座标系统，再由x,y座标系统，转换与网格间系，为了便于以后输出资料的分析，故应确实掌握网格之位置及其相对应的x,y 座标。建议在敏感区域使用较密之网格，其它地 区则使用较疏之网格，刚开始跑程式时， 不宜使用网格太大的网格数目，因尽量使 网格总数少于1000，以节省时间。4.2 FLAC输入程式编写步骤FLAC程式编写顺序需依实际工程进行之逻 辑步骤建立。如欲分析开挖题目，则需先求出未挖前之应力分布，再以此应力分布求开挖后之土壤变形及是否会崩毁，如因开挖后造成崩毁， 则FLAC程式将无法继续，萤幕将出现ERROR 讯息，如 Bad Geometry Zone ，-FLAC程式前几行之顺序为o Config_ o Grid_ o Model _ 以后各行编写之顺序，则建议以下之步骤：1. 求起始之应力平衡 (1)建立x,y座标与网格之关系，建议使用Gen指示：Gen x1,y1 x2,y2 ,x3,y3 x4,y4 i=i0,i1 j=j0,j1 详细指令参见使用手册，FLAC程式可自动产生x,y 座标与网格之关系，但由于产生之网格座标不易控 制，将对其它之工作产生负面影响，故依使用前步骤 4所建立之关系，将网格依其疏密程度需要之不同， 实际控制网格之座标。 (2)设定材料性质：prop (3)设定外力：Set Grav, Apply Pressure, ini sxx, Syy(4)设定边界条件：fix, free(5)求起始之应力平衡：solve (6)储存：Save 求工程之影响 求出区域内之应力分布情况后，再依工程之流程及步骤阶段执行各工程进行过程之影响，建议使用以下之步 骤： (1)叫出起初之应力平衡：re_ .sav (2)设定新的材料性质：model,prop (3)设定新的支撑性质：struct(4)设定新的外力 (5)设定边界条件 (6)求工程时之应力平衡 (7)储存 五、分析结果之印出及绘制FLAC V.3.4 在分析成果之绘制上，较先前之各版本 有一明显及方便的设计，因为V.3.4（CONSOLE）版 本。5.1 分析成果绘图a.直接绘图（不存图档）flac：set plot windowsflac：plot pen GR YD(欲画出格网及y向变位)说明:输入上述两行指令后，连结之印表机会直接印出图形， 约占半页A4之纸张，如欲印出全页，则须设定印表机横向 列印。 b.绘图(存图档)flac：set out YD.EMF(设以下要画之内容档名为.emf)flac：set plot emf color(设彩色印制)flac：plot pen GR YD(画出之内容为格网及Y向变位) 说明：YD.EMF可用Word软体叫出并绘图 5.2 印出分析内容或成果数据flac：set log YD.TXT(YD.TXT为所要列印内容之档名)flac：print yd i=1,10 j=1,5(设要印出y向位移量)flac：set log off说明：输入上述三个指令后，可用Word，Nodepad等软体叫出YD.TXT并列印。六、分析范例边坡稳定分析本范例共分析三种情况即 (1)粒性土壤 C=0 ；s13.sav (2)C土壤C0 0；s14.sav (3)考虑水位线 ；s15.sav 其分析网格之建立如下二图所示 01 title02 SLOPE UNDER GRAVITATIONAL LOAD03 grid 20,10 04 ;Mohr-Coulomb model 05 m m06 ;soil properties-note large cohesion to force initial elastic 07 ;behavior for determining initial stress sate. This will prevent08 ;slope failure when initializing the gravity stresses 09 prop s=.3e8 b=1e8 d=1500 fri=20 coh=1e10 ten=1e1010 ;warp grid to form a slope: 11 gen 0,0 0,3 20,3 20,0 j 1 4 12 gen same 9,10 20,10 same i 6 21 j 4 1113 mark i=1,6 j=4 14 mark i=6, j=4,1115 model null region 1,1016 ;displacement boundary conditions 17 fix x i=1 18 fix x i=21 19 fix x y j=120 ;apply gravity 21 set grav=9.8122 ;displacement history of slope23 his ydis i=10 j=1024 ;solve for initial gravity stresses 25 slove26 ;save initial stae27 save sll.sav28 ;reset displacement components to zero29 ini xdis=0 ydis=030 ;set cohesion to 0 31 prop coh=0 32 ;use large strain logic33 set large34 step 20035 plot hold bo dis xvel min-2.25e-4 max 0 int 2.5e-5 zero 36 save s12.sav37 step 80038 plot hold bo dis xvel min 1e-3 max 0 int 2e-4 zero39 save s13.sav40 rest sll.sav41 ini xdis=0 ydis=042 prop coh=le4 tens 0.0 43 set large 44 solve 45 save sl4.sav46 ini xdis=0.0 ydis=0.0 47 ;install phreatic surface in slope 48 water table 1 den 100049 table 1 (0,5) (6.11,5) (20,9) 50 def wet den 51 loop i (1,izones) 52loop j (1,jzones) 53if mode(i,j)1then54xa=(x(i,j)+x(i+1,j)+(i+1,j+1)+x(i,j+1) 55xc=0.25*xa 56ya=(y(i,j)+y(i+1,j)+y(i+1,j+1)+y(i,j+1) 57 yc=0.25*ya58 if yctable (1, xc) then59density(i,j)=180060end if61end if62 end loop 63 end loop 64 end 65 wet den 66 apply press 2e4 var 0 2e4 from 1,4 to 6,667 plot den blo water apply68 step 6000 69 sclin 1 19 0 19 10 70 plot hold bou vel pp 71 save s15.sav72 return 程式说明 Line 1 Line 2 本分析之抬头Line 3 建立20 10之格网 Line 4 以；为开始之指令为说明内容，FLAC不读；以后之指令 Line 5 分析模式为模耳库伦，即Model MohrLine 6 Line 8 说明：底下之土壤系数中之C值及张力强度值故意放大， 以免土壤在初始状况时即产生破坏Line 9 土壤性质 properties Shear Modulus G=0.3 108 N/m2 Bulk Modulus B=1 108 N/m2 density =1500 kg/m3 fri friction angle =20。Coh Cohesion C=1 108 N/m2 Ten Tension Strength T=1 108 N/m2G=E/2(1+v) k=E/(3(1-2v)Line 10 Line 12格网座标化 Line 13 Line 14以mark指令将格网区域化为两部份Line 15挖除由mark区分后包含Region(1,10)的那部份将边坡 制造出来Line 16 Line 19 边界条件之设定 i=1 及 21 之边界只容许上下位移，最 底部之边界则固定 Line 20 Line 21 加入重力于整个网格系统Line 22 Line 23 选择格点(10,10)检核其y向之位移历线Line 24 Line 25 开始计算求取答案Line 26 Line 27 储存初始平衡之条件。唯一之力为重力。此时之土体内 应力为边坡既有之情况Line 28 Line 29 将x向及y向之位移归零。因边坡在未被碰触之前其位 移应为零，但是其应力为力入重力后之自然现象，不应 归零 Line 30 Line 31 探讨CASE 1，纯粹性土壤之情况。将土壤之property 中 之Cohesion 设为零，其余参数不变 Line 32 Line 33 使用大应变之模式，则每一个step其格网座标自 动更新Line 34求解。除了使用slove之外，可令step=n，Line 35 绘出位移的等高线Line 37 设求解之step=800，以便观察800个step后之情况。 因本CASE中土壤之C=0 =20，很明显此边坡会破坏， step=800已够发现此现象。Line 40再叫出 sll.sav 即初始平衡况以便进行CASE 2 之分析Line 41再设X及Y向位移为零 Line 42 变更Cohesion为1 104 N/m2Line 47 Line 49 设定水位线为table 1，水的密度是1000 kg/cm3 Line 50 Line 64 使用FLAC去定义function wet_ den，此自设之 function乃在定义水位线下之单位重为 1800/cm3 。 注意：有几个 loop 就要有几个 end loop，有几个if 就要有几个end if，结尾必定要以end结束。 Line 65 前面之 define wet_ den 是在定义FLAC function之内 容。而此行才是要命令FLAC在此时执行此function .必需先定义出function之内容才能令FLAC执行。 Line 66 加外力于边界上，如在边界内部则使用interio指令。 此行即沿著i=1,4在j=6之边界加上2e4 Newton渐消 减为0之外应力。 Line 69定出画图之扫瞄线scan line，此行定出一条扫瞄线 连接(19,0)(19,10)两点，只要图中之等高线与此扫瞄 线相交，即会标出，方便看图时之参考。Line 72 跳出FLAC执行档 新关于flac的几点说明自从FLAC3D 随着时间的发展模拟了非线性体系，与常规的有限元程序在计算结束时生成结果相比，它的结果的解释说明更难了。这有一些指示器可以用来评估数值模型的状态例如，不管体系是稳定的，不稳定的，还是处于稳定的流变状态。各种各样的指示器用法如下所述。9.1不平衡力每个网格顶点最多由八个区域包围，这些区域对网格顶点施加力。在平衡状态，这些力的代数和几乎为零（也就是说，网格顶点一边的力几乎与另一边的力平衡）。如果不平衡力接近一个非零恒定值，那么这表示模型失败或进入了流变状态。在计算过程中，最大不平衡力由所有的网格决定，在显示屏上可以不断看到这个理。也可以把它保存为一个记录，并在图表里看到。对于评估模型的状态，不平衡力是很重要的，但是它的量级要与网格内典型的内力量级作比较。也就是说，有必要知道是什么产生了这个“小”力。网格顶点的典型内力可以根据增加与力垂直方向的应力来找到，要取网格重要区域内的典型值。用最大不平衡力与典型内力的比值表示R，表示为百分数，R的值从不会减小到零。但是，根据要求的精确度，1%或0.1%都可能被认为是达到了平衡（例如，在顺序执行的中间阶段可能R = 1%是足够好的，但是在做报告或论文时最终应力或位移分配可能要用R = 0.1% ）。注意R值小只表示所有网格顶点的力都平衡。但是可能会发生没有加速度的流变。为了区分流变和“真实”平衡，可能要检查另外的指示器，如下面所述。9.2网格顶点的速度评定网格的速度有两种方法，一是绘制整个区域的速度图（用PLOT vel 命令），一是选择网格中的一些关键点并记录它们的速度(HIS gp xvel, yvel或zvel)。两种绘图方法都有用。在最后阶段如果速度记录显示为水平线，那么表明达到了稳定状态。如果它们都收敛于趋近零（于它们的初始值相比），那么已经达到了绝对的平衡；如果有记录收敛于趋近零的值，那么与记录相应的网格顶点进入了流动状态。如果一个或多个速度记录图显示上下波动，那么系统可能出现了瞬时现象。注意速度由许多位移单元除以时步表示。速度向量图却很难说明，因为速度的大小和图样都很重要。只要网格顶点受力，速度就决不会减小到零。速度的大小应该与通过执行很多时步（如1000步）生成的位移有关。举例，如果体系的当前位移是1 cm ，且速度图中的最大速度是10-8 m/时步 ，那么1000时步会产生一个10-5 m, 或 10-3 cm 的位移增量，即当前位移的0.1%。这种情况下，即使速度在一个方向看起来是“平滑的”，那么也可以说体系是平衡的。更多的情况是，速度向量的方向很随意，（或几乎随意），大小可能很随意，（或几乎随意）。当网格顶点力的改变小于计算机的精度小数点后六位时，会发生上面那种情况。一个低振幅的速度场是没有流变的平衡状态的可靠的指示器。如果速度场的向量是一致的（也就是说，有一个系统的样式），且大小很大（使用前面所述的准则），那么体系家发生了流动或还在发生弹性变化（例如，发生了弹性阻尼振动）。为了确认流变是否在继续，应该检查一个塑性指示器图，如下所述。不管怎样，如果有弹性振动，那么就应观察速度大小，以显示这些运动是不是重要。可能会看到表面上看起来有意义的模式，但如果振幅比较低，那么运动没有物理意义。9.3塑性指示器对于FLAC3D 中的塑性模型，可用PLOT block state 命令显示塑