FLAC-FLAC3D基础与应用 ppt课件

GeoHohaiFLAC/FLAC3D基础与应用河海大学，2011年11月24日*ppt课件ppt课件*GeoHohai几个问题 什么是FLAC？ 为什么要用FLAC？ FLAC能做什么？ FLAC为何这么流行？ 怎么学FLAC？*ppt课件ppt课件GeoHohai什么是FLAC？ FastLagrangianAnalysisofContinua*ppt课件ppt课件GeoHohai为什么要用FLAC？Abaqus：745条Plaxis：80条ADINA：310条*ppt课件ppt课件GeoHohaiFLAC能做什么？ 岩土工程中的绝大多数问题 土力学、岩石力学、防灾减灾、隧道、地下空间等 采矿工程中的大部分问题 水工结构中的部分问题 结构工程国际通用的岩土工程专业分析程序*ppt课件ppt课件GeoHohaiFLAC为何这么流行？ CharlesFairhurst美国工程院、瑞典皇家工程院院士，国际岩石力学学科和岩石力学学会创始人之一，历任国际岩石力学学会主席和副主席，国际岩石力学学会Muller奖、美国岩石力学学会终生成就奖获得者。 PeterCundall美国工程院、英国皇家工程院院士，国际资深计算岩石力学学家。 *ppt课件ppt课件GeoHohai怎么学习FLAC？问：好学不？答：好学！问：难不？答：难！*ppt课件ppt课件GeoHohai安排 第一讲：基本介绍、静力分析、前后处理 第二讲：接触面、FISH语言、流固耦合分析 第三讲：动力分析、自定义本构、结构单元 第四讲：FLAC（2D）基本介绍与应用实例 讨论*ppt课件ppt课件GeoHohai第一讲FLAC3D基本介绍、静力分析、前后处理*ppt课件ppt课件GeoHohai软件介绍 FastLagrangianAnalysisofContinua 美国Itasca咨询公司开发2D程序(1986) 1990年代初引入中国 有限差分法(FDM) 3D版本：DOS版→2.0→2.1→3.0→3.14.0 2D版本：DOS版4.05.06.07.0*ppt课件ppt课件GeoHohaiFLAC 大应变、小应变计算模式. 丰富的本构模型、提供自定义的本构模型功能 接触面可以模拟不同材料的接触 流固耦合实现土体的固结与渗流 拥有各种功能的结构单元类型，模拟土与结构的相互作用 强大的动力分析功能. 流变分析，拥有粘弹性模型和粘塑性模型 热力学分析.Shearstrainratecontours*ppt课件ppt课件*GeoHohaiFLAC3D 与FLAC类似，是FLAC的三维版本 与FLAC拥有相同的优点upstreamdownstream*ppt课件ppt课件*GeoHohai基本特点内置材料模型连续介质非线性，大应变模拟显式解题，为不稳定物理过程提供稳定解界面或滑动面用来模拟可产生滑动或分离的离散面，从而模拟断层，节理或摩擦边界内置材料模型丰富：零模型,三个弹性模型(各向同性，横观各向同性和正交各向异性),八个朔性模型(德鲁克-布拉格,摩尔-库伦,应变硬化/软化，单一节理，双线性应变硬化/软化单一节理,双屈服，修正剑桥粘土，霍克-布朗)*ppt课件ppt课件GeoHohai可选模块可选模块包括:热力学,热-力学耦合,热-流体-力学耦合包括热传导和对流;粘弹,粘朔性(蠕变)材料模型;动力学分析,并可以模拟静边界和自由域使用C++定义自己的模型核废料储存中的热力学研究问题*ppt课件ppt课件GeoHohaiFLAC/FLAC3D基本原理 FLAC/FLAC3D利用有限差分，显示方案，动态松弛方法模拟连续体的非线性力学行为：即使对准静态问题，程序仍然求解完整的动力学方程。这种方法的好处在于可以为物理非稳定过程例如塌方提供稳定解；在“松弛”方法中，使用阻尼来吸收动能以模拟系统的“静态”反应。这种方法可以用比其它方案如解矩阵法更为真实有效地模拟塌方问题。*ppt课件ppt课件GeoHohaiLagrangian法 源自流体力学中的拉格朗日法 跟踪流体质点的运动状态 跟踪固体力学中结点，按时步用Lagrangian法研究网格节点的运动 节点和单元随材料移动，边界和接触面与单元的边缘一致 固体力学大变形理论法国数学家、物理学家拉格朗日*ppt课件ppt课件GeoHohai混和离散技术FLAC混和离散+/2=每个为常应力/应变:体积应变由整个四边形算出.应变偏量则有两个三角形和分别算出(混合离散过程)解题过程中网格坐标按照“拉格朗日方式更新”(网格随材料移动),且为显式(一个时步内局部变化不会影响邻域)*ppt课件ppt课件GeoHohaiFLAC3D混和离散+/2=*ppt课件ppt课件GeoHohaiFLAC3D混和离散 结构域离散为可由四面体单元组合形成的五面体或六面体等单元; 以为基本单元(常应力、常应变);体应变的计算：;偏应变的计算：.*ppt课件ppt课件GeoHohai动态松弛动态松弛法 在动态松弛法中，网格点根据牛顿运动定律运动.网格点的速度与该点的不平衡力呈正比.这种求解方法所决定的一系列位移将把系统带入平衡状态，或表明破坏模式. 在动态松弛法中有两个因素很重要: 时步的选择 阻尼效应*ppt课件ppt课件GeoHohai显式算法显式解与隐式解的比较显式,逐时推进隐式,静态1.无需进行反复迭代来实现非线性本构关系.2.类似问题求解时间呈N3/2规律增长3.物理非稳定性不会引起数值不稳定性.4.因为无需储存矩阵，用较小内存即可模拟大尺度问题.5.对大位移、大应变问题同样适合，无需额外的计算.1.需进行反复迭代来实现非线性本构关系2.类似问题求解时间呈N2甚至N3规律增长.3.难以模拟物理非稳定性问题.4.需存储刚度矩阵，需克服相关的带宽问题，需要的内存较大.5.对大位移、大应变问题需进行大量的计算.*ppt课件ppt课件GeoHohaiNewFeaturesinFLACVersion6.0 使用IntelFortrancompiler拥有更快的计算速度 自动网格重画功能，解决bad-geometry问题. 新的模拟颗粒土材料的硬化模型 更新的通用网格生成工具*ppt课件ppt课件*GeoHohaiNewFeaturesinFLAC3DVersion3.1 多处理器的并行计算功能 新结构单元类型“EmbeddedLiner”提供两个方向的接触作用，可以很好地模拟挡土墙 对四面体单元采用新的混合离散方法“NodalMixedDiscretization”提供塑性问题更精确的解答 64位程序 包含命令手册、FISH手册和应用实例的帮助*ppt课件ppt课件GeoHohaiNewFeaturesinFLAC3DVersion4.0 模拟颗粒状材料的硬化模型 自动网格重画功能，解决bad-geometry问题. 改进的interface 更快的渗流计算 更新的动力计算功能*ppt课件ppt课件GeoHohaiLagrangian格式动量平衡方程F(t)m牛顿运动定律对于连续体在静力平衡条件下，加速度项为0，方程变为平衡方程*ppt课件ppt课件GeoHohai自由落体的模拟G=mgS=1/2gt2=20m命令流：configdyngenzonbrisize111inixmul0.1ym0.1zm0.1modelelaspropbulk3e8shear1e8inidens1000setgrav00-10solveage2*ppt课件ppt课件GeoHohai自由落体的模拟(movie)*ppt课件ppt课件GeoHohaiFLAC3D中模型术语节点gridpoint：节点zone：单元boundary：边界*ppt课件ppt课件GeoHohaiFLAC3D的求解过程*ppt课件ppt课件GeoHohaiFLAC3D中的本构模型 开挖模型null 3个弹性模型 各向同性弹性 横观各向同性弹性 正交各向同性弹性 8个弹塑性模型 Drucker-Prager模型、Morh-Coulomb模型、应变硬化/软化模型、遍布节理模型、双线性应变硬化/软化遍布节理模型、修正剑桥模型和胡克布朗模型*ppt课件ppt课件GeoHohaiFLAC3D中的本构模型*ppt课件ppt课件GeoHohai一个最简单的例子genzonbrisize333;建立网格(前处理)modelelas;材料参数propbulk3e6shear1e6inidens2000;初始条件fixzranz-.1.1;边界条件fixxranx-.1.1fixxranx2.93.1fixyrany-.1.1fixyrany2.93.1setgrav00-10solve;求解appnstr-10e4ranz3x12y12solveploconzd;后处理 切片功能RUNFLAC3D*ppt课件ppt课件GeoHohai分析问题的过程建立网格初始条件边界条件初始应力平衡外荷载求解前处理后处理*ppt课件ppt课件GeoHohaiFLAC3D的文件格式保存文件(*.sav)–含有所有状态变量和用户定义条件的二进制文件数据文件(*.dat)–数据文件由用户创建的一种ASCⅡ格式的文件，它包括一系列的用于描述所分析问题的FLAC3D命令FISH文件(*.fis)–FISH程序文件FLAC3D文件(*.flac3d)–FLAC3D的网格信息文件历史文件(*.his)–记录输入输出历史值的文件图形文件–图形文件(各种标准格式)电影文件(*.dcx)–AVI或PCX图像文件，这些图像文件可以当作电影放映*ppt课件ppt课件GeoHohai初始应力的生成 为什么要单独列出？ 分析过程中出现的很多问题都与初始应力是否合理有关 手册中的例子五花八门 是所有后续分析的基础！！ 生成方法 弹性求解 更改强度参数的弹塑性求解 设置初始应力的弹塑性求解 存在水压力的初始应力生成 水下建筑的初始应力生成*ppt课件ppt课件GeoHohai弹性求解genzonbrisize112melaspropbulk3e7shear1e7fixzranz0fixxranx0fixxranx1fixyrany0fixyrany1inidens2000setgrav00-10solveStep=162sz=-40e3sx=-21.54e3*ppt课件 模型尺寸 单元数量 密度 K G u 1×1×2(m3) 1×1×2 2000 30MPa 10MPa 0.35ppt课件GeoHohai更改强度参数的弹塑性求解genzonbrisize112modelmohrpropbulk3e7shear1e7c1e10f15ten1e10fixzranz0fixxranx0fixxranx1fixyrany0fixyrany1inidens2000setgrav00-10solvepropbulk3e7shear1e7c10e3f15ten0solveStep=163sz=-40e3sx=-21.54e3Or:solveelastic*ppt课件 模型尺寸 单元数量 密度 K G c f u 1×1×2(m3) 1×1×2 2000 30MPa 10MPa 10kPa 15 0.35ppt课件*GeoHohai设置初始应力的弹塑性求解genzonbrisize112modelmohrpropbulk3e7shear1e7c10e3f15ten0fixzranz0fixxranx0fixxranx1fixyrany0fixyrany1inidens2000iniszz-40e3grad0020e3ranz02inisyy-20e3grad0010e3ranz02inisxx-20e3grad0010e3ranz02setgrav00-10solveStep=0sz=-40e3sx=-20e3*ppt课件 模型尺寸 单元数量 密度 K G c f u 1×1×2(m3) 1×1×2 2000 30MPa 10MPa 10kPa 15 0.35ppt课件GeoHohai存在水压力的初始应力生成(1)genzonebricksize112modelmohrpropbulk3e7shear1e7coh10e3fri15ten0fixzranz0fixxranx0fixxranx1fixyrany0fixyrany1inidens2000ranz01inidens1500ranz12iniszz-35e3grad0020e3ranz01inisyy-22.5e3grad0015e3ranz01inisxx-22.5e3grad0015e3ranz01iniszz-30e3grad0015e3ranz12inisyy-15e3grad007.5e3ranz12inisxx-15e3grad007.5e3ranz12inipp10e3grad00-10e3ranz01setgrav00-10solveStep=0sz=-35e3sx=-22.5e3rd=rs–n×s×rf*ppt课件 模型尺寸 单元数量 饱和密度 K G c f u 水位线 孔隙率 1×1×2(m3) 1×1×2 2000 30MPa 10MPa 10kPa 15 0.35 1m 0.5ppt课件GeoHohai存在水压力的初始应力生成(2)configfluidgenzonbrisize112modelelaspropbu3e7sh1e7inidens1500modelfl_isoinifdens=1000fmod0proppor0.5setgrav00-10waterdens1500watertableface001,011,111,101inipp10e3grad00-10e3ranz01iniszz-30e3grad0015e3ranz12...setfluidoffsolveStep=142sz=-40e3sx=-24e3单元数较少产生的误差*ppt课件 模型尺寸 单元数量 饱和密度 K G c f u 水位线 孔隙率 1×1×2(m3) 1×1×2 2000 30MPa 10MPa 10kPa 15 0.35 1m 0.5ppt课件GeoHohai水下建筑的初始应力生成genzonbrisize112modelmpropbulk3e7shear1e7c10e10f15ten1e10fixzranz0fixxranx0fixxranx1fixyrany0fixyrany1inidens2000ranz02iniszz-50e3grad0020e3ranz02inisyy-40e3grad0015e3ranz02inisxx-40e3grad0015e3ranz02inipp30e3grad00-10e3ranz02appnstress-10e3ranz2setgrav00-10solveStep=0sz=-50e3sx=-40e3*ppt课件 模型尺寸 单元数量 饱和密度 K G c f u 水位线 1×1×2(m3) 1×1×2 2000 30MPa 10MPa 10kPa 15 0.35 3mppt课件GeoHohai前后处理*ppt课件ppt课件GeoHohai前后处理 基本前后处理 命令操作 菜单操作 dd&dip attach&merge 外界模型的导入 复杂模型的网格检查*ppt课件ppt课件GeoHohaiFLAC3D的前后处理 命令驱动(推荐) 程序控制 图形界面接口 计算模型输出 指定本构模型及参数 指定初始条件及边界条件，指定结构单元 指定接触面 指定自定义变量及函数(FISH) 求解过程的变量跟踪 进行求解 模型输出*ppt课件ppt课件GeoHohai菜单驱动(计算模式)命令栏*ppt课件ppt课件GeoHohai菜单驱动(Plot)*ppt课件ppt课件GeoHohaiFLAC3D的前处理FLAC3D网格生成的关键特征： FLAC3D是命令驱动. 使用FLAC3D内置基元进行形状组合可形成复杂网格. 用户自定义FISH函数可以用来修改基元网格以创建更为复杂的网格. 第三方软件导入.*ppt课件ppt课件GeoHohaiFLAC3D网格基元*ppt课件ppt课件GeoHohai*ppt课件ppt课件GeoHohai块体径向柱体柱状交叉*ppt课件ppt课件GeoHohaiTecplot的后处理*ppt课件ppt课件GeoHohai*ppt课件ppt课件GeoHohaiTecplot后处理*ppt课件ppt课件GeoHohai*ppt课件ppt课件GeoHohai*ppt课件ppt课件GeoHohai*ppt课件ppt课件GeoHohaiTecplot的云图和等值线Tecplot中云图的效果与flac3d的shadeon的效果差不多,但是其出三维等值线的功能是卓越的.*ppt课件ppt课件GeoHohai*ppt课件ppt课件GeoHohaiTecplot的切片功能Tecplot的slice切片功能.和CAD一样,可以任意切剖面出图.最大优点是可以几个剖面同时出图.*ppt课件ppt课件GeoHohai整体和截面的数据以单元形式存在列表中,可自由选择和组合多个单元出图.*ppt课件ppt课件GeoHohaidd和dip Simwe.Com上有近百的讨论贴 地质上的概念，倾向和倾角 建议用ori和norm代替 后处理的切片功能 plosetplaneori(*,*,*)norm(*,*,*) ploconzdisplane 概念清晰，使用快捷 其他应用 reflect网格 建立水面xyzorinorm*ppt课件ppt课件GeoHohaiattach和merge attach 用于连接次节点(sub-grid) 允许网格存在一定随意性 要求成比例(1:2,1:3等) 可用于模型的检查 慎用 merge 用于节点的合并 外来模型导入的精度差异 tolerance的设置sub-grid4:25:2tol121merg*ppt课件ppt课件GeoHohai外界模型的导入 复杂网格的生成难度大 接口编写 不同软件之间的精度差异 Group的定义 采用.flac3d的文件格式 .flac3d文件的格式 G11.0e+001.0e+001.0e+00 ZB812345678 ZGROUPSoil123 impgrid&expgrid 仅限于网格3wzones3sec*ppt课件ppt课件GeoHohai复杂网格的检查 FLAC3D生成的复杂网格 attachface 无接触面时检查整体模型是否存在sub-grid 有接触面时给定范围进行检查 其他软件导入的复杂模型 网格划分的检查 genmerge 弹性模型modelelastic求解 “独立”节点 “畸形”单元*ppt课件ppt课件GeoHohai模型的检查 FLAC本身的Check功能十分有限 错误提示很少 十分开放的工作平台(自由落体) 检查的基本步骤 网格检查(如前所述) 边界条件检查 速度约束条件 plogpfixredsk 荷载条件 plofapredsk 模型检查 模型赋值 ploblockmodel 参数赋值 ploblockprop*** 初始应力检查(如后所述)*ppt课件ppt课件GeoHohai初始应力检查办法 初始应力的计算时间不会“非常长” 经常检查模型的响应 ploconszz(syy,sxx) 应力场 ploconzdis(ydis,xdis) 位移场 ploblosta 屈服状态 plogpfixredsk 速度约束条件 plofapredsk 体力 plohist(unbal) 不平衡力 plointerfacenstress(sstress) 接触面单元*ppt课件ppt课件GeoHohai第二讲FLAC3D接触面、FISH语言、流固耦合分析*ppt课件ppt课件GeoHohai接触面（interface）*ppt课件ppt课件GeoHohai接触面单元 应用范围 原理 建模方法 推荐方法 复杂内部接触面的设置方法 参数选择 单桩承载力分析 挡土墙的接触面设置 思路：未知问题的分析方法*ppt课件ppt课件GeoHohai接触面单元的用途 网格不连续 岩体介质中的解理、断层、岩层面 地基与土体的接触 箱、槽及其内充填物的接触 空间中无变形的固定“障碍”*ppt课件ppt课件GeoHohai接触面的原理 三角形单元(无厚度!) 参数较多 三种工作模式 粘结界面 粘接滑移 库伦滑动*ppt课件ppt课件GeoHohai接触单元模型的建立(1) 关键要形成同一位置的两个节点(面) “移来移去”(推荐) 建两个分开的模型 建立接触单元 通过INI*add使模型接触 注意dist的含义 NOmerge,NOattach!!!接触面dist1234*ppt课件ppt课件GeoHohai接触单元模型的建立(2) “导来导去” 利用expgrid,impgrid命令进行网格导出与导入 配合DELETE命令 适于内部接触面的建立，或 其他前处理工具建立的网格*ppt课件ppt课件GeoHohai接触单元模型的建立(3) GENseparate INTERFACEwrap 指定正确的group*ppt课件ppt课件GeoHohai接触面参数的确定 虚构的为了合并节点而设置的接触面 Kn=ks=10*max[(K+4/3G)/Dzmin] 真实的刚性接触面 如料仓下料 c,D,Tension重要，kn,ks不重要 真实的柔性接触面 断层；水力劈裂材料 试验得到参数 对于kn,ks：岩石断层10~100MPa/m(粘土);100GPa(岩石) 反分析方法：通过断层中岩石的变形与原岩的变形*ppt课件ppt课件GeoHohai单桩承载力分析 软土地基 bulk1.6878E6 shear3.6167E5 coh15E3 fric12 dens1.73E3 桩体 bulk5e9 shear3.75e9 dens2.5e30.5m8m10m20m*ppt课件ppt课件GeoHohai计算过程施加桩顶荷载计算结果*ppt课件ppt课件GeoHohai影响因素7.50E+03fricgk3coh/0.70.7fric10g10k23.00E+0420100e100e1cohfrickskn水平因素*ppt课件GeoHohai计算工况设计1233(9)3123(8)2313(7)2132(6)1322(5)3212(4)3331(3)2221(2)1111(1)cohfrickskn水平因素工况*ppt课件GeoHohai计算结果Ks取1Ks取2Ks取396.9%1233(9)22.4%3123(8)54.1%2313(7)96.9%2132(6)7.1%1322(5)37.8%3212(4)96.9%3331(3)42.9%2221(2)37.8%1111(1)Ra差异度cohfrickskn水平因素*ppt课件GeoHohai最优方案2222111122101122最优方案17.3%6.8%72.8%11.9%极差52.4%52.7%96.9%57.8%k364.6%59.2%24.1%47.3%k247.3%52.4%43.2%59.2%k1cohfrickskn水平因素*ppt课件GeoHohai合理步骤单桩分析简单网格接触面参数多次试算理论、实测加密网格接触参数理想结果群桩分析Pile结构单元单元参数理想结果Pile结构单元单元参数理想结果*ppt课件ppt课件GeoHohai挡土墙的接触面设置 对于未知问题的分析思路 3个独立的接触面 3个同ID的接触面 2个独立的接触面 2个独立的接触面并进行底部merge挡墙土体123WallSoil*ppt课件ppt课件GeoHohai挡土墙的接触面设置(2)Z=5.16cmZ=1.14cmZ=∞不能初始平衡Z=0.35cm？*ppt课件ppt课件GeoHohai挡土墙的接触面设置(3)X=1.95cmX=1.91cmX=1.47cm？*ppt课件ppt课件GeoHohai挡土墙的接触面设置(3) 前两种方法的差别实质 3个ID的独立接触面在相同位置产生互不影响的两个节点 共同ID的接触面在相同位置自动设置为1个节点 最终的结论需要您自己去判断！3interfaces,2IDs1interface,1ID*ppt课件ppt课件GeoHohaiFISH语言*ppt课件ppt课件GeoHohaiFISH语言简介 软件自带的编程语言 是否一定要学？ 视情况而定，需要时查询FISH变量即可 语法简单xxx……end_xxx 注意事项 与FLAC本身的关键字冲突 保留字不可缩写 变量可不定义，因此注意检查程序 printfish table,extra等命令使用*ppt课件ppt课件GeoHohai一个最简单的FISH程序defabcabc=1+2*3abcd=1.0/2.0endabcprintfish数据格式*ppt课件ppt课件GeoHohai函数与变量 都可以在FISH函数中进行赋值，赋值操作与常规的编程语言类似，按照运算符的优先级先后顺序来执行。 函数和变量的赋值遵守数据类型的，即整型的计算结果为整型，浮点型的计算结果为浮点型，因此读者在进行除法运算、开方运算时都需要将数据类型设置为浮点型，数字尽量使用小数点以保证运算正确。 变量和函数名的命名规则是不能以数字开头，不能含有中文，并且不能包含如下的字符。.,*/+-ˆ=<>#()[]@;’" 变量和函数名不能与FLAC3D、FISH的保留字相冲突，不要采用过于简单的单词，比如a，hist等，这些都与保留字相冲突。 即使程序中存在与保留字相冲突的变量，FLAC3D也不会提供任何提示，所以提醒读者在编制FISH程序时尽量使用较长的、复杂的变量和函数名。*ppt课件ppt课件GeoHohai函数与变量 对变量进行赋值时，不能使用当前函数的函数名放在“=”的右边，比如采用下面的定义abcd=abc+1.0在FISH程序执行时会提出错误，因为这样会形成递归调用，这种调用方式在FISH程序中是不允许的。 变量和函数的作用是全局的，在命令中的任何地方修改变量的值都会立即生效，因此在实际应用中尽量避免不同的函数中含有相同的变量，因为这样可能会造成赋值错误，并难以检查。 在FLAC3D中可以用如下的命令来引用FISH函数和变量 PRINT用于查看函数和变量的数值； HISTORY命令可以对函数和变量的数值进行记录； SET命令用于变量的赋值。*ppt课件ppt课件GeoHohai主要语句 选择语句CASEOF表达式…默认语句CASEn1…表达式的值为n1时的语句CASEn2…表达式的值为n2时的语句ENDCASE*ppt课件ppt课件GeoHohai主要语句 条件语句IF条件表达式[THEN]…[ELSE]…ENDIFFISH中条件运算符没有“并”、“或”、“否”这样的符号ifaa>1.0ifaa<2.0执行语句endifendif表达“1<aa<2”的条件*ppt课件ppt课件GeoHohai主要语句 循环语句LOOPvar(exp1,exp2)…ENDLOOP或者LOOPWHILE条件表达式…ENDLOOPp_z=zone_headloopwhilep_z#null;语句p_z=z_next(p_z)endloopp_gp=gp_headloopwhilep_gp#null;语句p_gp=gp_next(p_z)endloop单元遍历语句节点遍历语句*ppt课件ppt课件GeoHohai主要语句 命令语句COMMAND…（FLAC3D命令）ENDCOMMAMD*ppt课件ppt课件GeoHohai应用实例： 土体的模量随小主应力变化defE_modifyp_z=zone_headd_k=704d_n=0.38d_pa=101325.0;//标准大气压loopwhilep_z#nullsigma_3=-1.0*z_sig1(p_z)E_new=d_k*d_pa*(sigma_3/d_pa)^d_nz_prop(p_z,'young')=E_newp_z=z_next(p_z)endloopendE_modify*ppt课件ppt课件GeoHohai应用实例 获得最大位移的大小及发生位置deffind_max_dispp_gp=gp_headmaxdisp_value=0.0maxdisp_gpid=0loopwhilep_gp#nulldisp_gp=sqrt(gp_xdisp(p_gp)^2+gp_ydisp(p_gp)^2+gp_zdisp(p_gp)^2)ifdisp_gp>maxdisp_valuemaxdisp_value=disp_gpmaxdisp_gpid=gp_id(p_gp)endifp_gp=gp_next(p_gp)endloopendfind_max_dispprintmaxdisp_valuemaxdisp_gpid*ppt课件ppt课件GeoHohaiFISH的编写习惯*ppt课件 第一步 第二步 defabcendabc defabcp_gp=gp_headloopwhilep_gp#nullp_gp=gp_next(p_gp)endloopendabc 第三步 第四步 defabcp_gp=gp_headloopwhilep_gp#nullcommandendcommandp_gp=gp_next(p_gp)endloopendabc defabcp_gp=gp_headloopwhilep_gp#nullcommandappnstress…endcommandp_gp=gp_next(p_gp)endloopendabcGeoHohaiFISH与建模两个圆形隧道的连接部分变直径的隧道部分*ppt课件ppt课件GeoHohaiFISH检查 主要采用PRINTfish命令 查看变量的赋值是否合理，主要检查值为0的函数和变量，因为FISH程序中一般定义的变量都有实际的意义，输出为0的变量很可能是与保留字相冲突的变量 （如a就是apply的保留字） 或者由于编写笔误产生的变量 （如数字0与大写字母O）。*ppt课件ppt课件GeoHohai流固耦合分析很难！*ppt课件ppt课件GeoHohai流-固耦合分析(单相流) 基本功能 理论框架 计算模式 渗流边界条件，初始条件 单渗流计算及渗流耦合计算*ppt课件ppt课件GeoHohai基本功能 渗流各向同性、各向异性 不同的渗流模型和属性 流体压力，涌入量，渗漏量和不渗水边界 抽水井、点源、体积源 饱和渗流可采用显式差分法、隐式差分法 非饱和渗流采用显式差分法 渗流-固体-热的耦合 流体和固体的耦合程度依赖于土体颗粒(骨架)的压缩程度，用Biot系数表示颗粒的可压缩程度。 循环荷载引起的动水压力变化和土体液化*ppt课件ppt课件GeoHohai理论框架 准静态Biot理论 多孔介质中遵循Darcy定律的单相渗流 描述多孔介质中流体渗流的变量 孔隙水压力，饱和度，特定排水向量的三个分量 质量守恒定律 达西定律 本构定律 考虑流体响应孔隙水压力改变，饱和度改变，体积应变改变和温度改变*ppt课件ppt课件GeoHohai有效应力计算 无渗流模式 不设置CONFIGFluid 孔压不改变 设置孔压分布 INITIALpp WATERtable WATERdensity WATERtableface SETgravity 手动设置干湿密度 渗流模式 设置CONFIGfluid 设置土体干密度 渗流模型 MODELfl_isotropic MODELfl_anisotropic MODELfl_null SETfluidoff setWATERbulk=0*ppt课件ppt课件GeoHohai渗流边界条件，初始条件 默认的边界条件是不透水边界 孔隙压力自由(不透水边界) 固定孔隙水压力(透水边界) 如：井 孔隙压力，孔隙率，饱和度和流体属性的初始分布可以用INITIAL命令或者PROPERTY命令定义。*ppt课件ppt课件GeoHohai单渗流计算及渗流耦合计算 时间比例 完全耦合分析方法 孔压固定分析(有效应力分析) 单渗流得到孔压分布 无渗流计算——孔压的力学响应 流-固耦合计算*ppt课件ppt课件GeoHohai单渗流得到孔压分布 用途：排水沟；抽水井；耦合计算 计算步骤 CONFIGfluidSETmechoff SETfluidimpliciton/off MODELfl_;PROP STEP;SOLVEage;SETfluidratio SETfluidoffmechon PROPbiot_c0(orINIfmod0)*ppt课件ppt课件GeoHohai无渗流计算——孔压的力学响应 不排水短期响应 两种分析方法：干法和湿法 干法：Ku=K+a2M 两种破坏形式 WATER或INI获得常孔压，不排水的c,φ(孔压改变较小) φ=0,c=cu(M>>K+4/3G) 湿法：耦合体系的短期行为 使用排水的K,c,φ 若SETfluidoff,Biot_mod(fmod)真实*ppt课件ppt课件GeoHohai力学过程的特征时间流体扩散过程的特征时间完全耦合分析方法 时间比例 短期行为(不排水) ts(分析时间)<<tc(耦合扩散时间) 忽略渗流影响 长期行为(排水)ts>>tc 施加扰动的属性 流体扰动：渗流可不与力学过程耦合 力学扰动：耦合等级取决于流固刚度比 流固刚度比*ppt课件GeoHohai流-固耦合计算 CONFIGfluid;M(Kf);K(渗透系数)真实，则FLAC3D默认耦合计算 Δp→Δev→s Δev→Δp 预估流/力特征时间 耦合计算前先达到一个平衡状态 SETfluidonmechoff;SETfluidoffmechon;STEP SETmechforce;SETmechsubstepnauto;SETfluidsubstepm(=1) STEP:渗流步足够小*ppt课件ppt课件GeoHohai流固耦合的计算方法 手动调整的STEP求解 SETfluidonmechoff STEP SETfluidoffmechon STEP 主从进程的SOLVE求解 SETmechforce SETmechsubstepnauto (从进程) SETmechsubstepm (主进程) SOLVEage 自动STEP求解 STEP*ppt课件ppt课件GeoHohai渗流问题(CONFIGfluid)分析步骤 时间比例(ts,tc) 稳态 不排水状态 相当 扰动类型 力学扰动 孔压扰动 流固刚度比Rk 是否>>>1 完全耦合模式 时间比例相当；力学扰动*ppt课件ppt课件*Ts分析的时间Tc扩散的时间GeoHohai心墙土坝的渗流(1)newconfigfluidsetfluidoffgenzonbrickp000-10size20110genzonbrickp0500p11500p2510p3905p41510p5915p61105p71115size1015groupsoilgroupdamranx57z-50groupdamranid201aid211aid221aid231aid241agroupdamranid202aid212aid222aid232aid242amepropbu3e7sh1e7inipp0grad00-10e3ranz0-10inidens2000modelfl_isoproppor0.5perm1e-10inifden1000ften-1e10inisat0.0ranz05modelfl_nullrangrodam;inipp0rangrodamfixzranz-10fixxranx0fixxranx20fixysetgrav10solvesaveelastic.sav网格模型初始孔压*ppt课件ppt课件GeoHohai心墙土坝的渗流(2)restelastic.savinixd0yd0zd0xv0yv0zv0appnstress-40e3grad0010e3ranz04x09solvesavepressure.sav竖向应力沉降*ppt课件ppt课件GeoHohai心墙土坝的渗流(3)restpressure.savsetfluidonmechoffinifmod2e3ften0.0rangrosoilinixd0yd0zd0xv0yv0zv0apppp40e3grad00-10e3ranz04x09apppp0ranz0x1520histid=10zoneppid215solveRatio=1Ratio=1E-5*ppt课件ppt课件GeoHohai荷载引起的地基土体的超孔隙水压力*ppt课件ppt课件GeoHohai计算文件;---applyloadslowly---deframpramp=min(1.0,float(step)/200.0)endapplynstress=-40e3histramprangex-.13.1z9.910.1;---fluidflowmodel---modelfl_isoinifmod2e9;---porepressurefixedatzeroatthesurface---fixpp0rangez9.910.1;---settings---setfloff;---test---step750[1][1]因为本例中没有设置初始应力，这里只进行了750步的求解。*ppt课件ppt课件GeoHohai计算结果*ppt课件ppt课件GeoHohai对主从进程法的讨论SETmechforce设置一个不平衡力的大小，达到这个不平衡力系统认为暂时达到平衡状态；SETmechsubstepnauto设置力学进程为从进程，在主进程每执行一步中必须执行n步，当系统达到平衡时也可以少于n步；SETfluidsubstepm设置流体进程为主进程。讨论1：对收敛准则进行对比分析不平衡力（force）：1E3、5E3、1E4、5E4不平衡力比（ratio）：1E-4、1E-3、1E-2、1E-1采用1E-3的收敛准则既可以满足计算流固耦合过程中的计算精度要求，同时又具有较高的计算效率*ppt课件ppt课件GeoHohai对主从进程法的讨论(2)讨论2：子步数的影响setmechsub100fluidsub10setmechsub10fluidsub10setmechsub1fluidsub1设置合理的子步数也很重要。设置过大，则会导致计算时间大大增加，过小又会造成计算结果的误差。*ppt课件ppt课件GeoHohai真空预压的简单模拟 孔压边界条件 ts>>tc 长期分析(排水) Rk>>1 骨架很软 孔压扰动 进行biot_mod调整砂层软土层粘土层PVD2m8m10mDatafile:*ppt课件GeoHohai数值分析过程(movie)*ppt课件ppt课件GeoHohaiFLAC3D非线性动力分析非常复杂！SaidbyProf.PeterCundall*ppt课件ppt课件GeoHohai为什么要用FLAC做动力分析？ FLAC可以模拟体系（土，岩石，结构，流体）受到的外部动力荷载（比如地震）或内部动力荷载（比如基础振动、爆炸）。 可以计算塑性引起的永久变形以及孔隙水压力的消散。 土动力学中常用的等效线性方法无法直接处理上述问题。*ppt课件ppt课件GeoHohai动力模拟的3个重要问题 动力荷载与边界条件 材料响应与阻尼 土体液化*ppt课件ppt课件GeoHohai动力荷载 动力输入的类型 加速度时程 速度时程 应力(压力)时程 力时程 APPLY INTERIOR(内部) TABLE FISH*ppt课件ppt课件GeoHohaiQuiet边界 静态(quiet,粘性)边界 LysmerandKuhlemeyer(1969) 模型边界法向和切向设置独立的阻尼器 性能 对于法向p波和s波能很好的吸收 对于倾斜入射的波和Rayleigh波也有所吸收，但存在反射 人工边界仍应当足够远*ppt课件ppt课件GeoHohaiQuiet边界应用 内部振动(如隧道中的列车振动问题) 动力荷载直接施加在节点上 使用Quiet边界减小人工边界上的反射 不需要FF边界 外部荷载的底部边界 软土地基上的地震荷载不适合用加速度或速度边界条件 使用应力条件t=-2Csrvs 地震底部输入的侧向边界 扭曲了入射波quietquietquiet*ppt课件ppt课件GeoHohaiFree-field边界 Cundalletal.(1980) 自由场网格与主体网格的耦合粘性阻尼器，自由场网格的不平衡力施加到主体网格边界上 设置条件 底部水平，重力方向为z向 侧面垂直，法向分别为x,y向 其他边界条件在APPLYff之前相当于一个阻尼器*ppt课件ppt课件GeoHohaiFree-field边界 APPLYff将边界上单元的属性、条件和变量全部转移ff单元上； 设置以后主体网格上的改动将不会被FF边界所响应 可存在任意的本构模型以及流体耦合(仅竖向) FF边界进行小变形计算，主体网格可大变形，FF边界上的变形要相对较小 存在attach的边界将不能设置FF边界 边界上的Interface将不能连续 动力边界设置需在FF边界设置之前*ppt课件ppt课件GeoHohaiFree-field边界与动力荷载 模型底部边界 fix——施加速度或加速度荷载——刚性边界 Free——施加应力时程荷载——柔性边界 对于软弱的地基不适合施加速度（加速度荷载），而应当施加应力荷载 Notethatthereisafactorof2becausetheinputenergydividesintoadownward-&upward-propagatingwave.*ppt课件ppt课件GeoHohai2.材料响应与阻尼 连续的非线性，表观模量随着应变的增大而降低 对所有循环应变等级均存在滞回特性，因此导致随着循环应变的增加阻尼比增大。阻尼是率相关的。 对于复杂波形的各个成分都产生阻尼。 剪切应变会产生的体积应变，相应的，随着剪应变循环次数的增加体积应变逐渐积累。*ppt课件ppt课件*GeoHohai材料响应土体在循环荷载作用下呈现出模量衰减和能量消散的特点，那么如何用非线性数值方法对其进行模拟呢？Nonlinearcharacteristicsofsoils(MartinandSeed,1979)*ppt课件ppt课件*GeoHohai试验得到的阻尼比、割线模量随循环剪应变的曲线rmalizedShearModulus,G/GmaxMid-RangeSandCurve(Seed&Idriss,1970)SandFillInland:Friction=32,hr=0.47,Go=440SandFillunderRockDike:Friction=30,hr=0.43,Go=440*ppt课件ppt课件*GeoHohai等效线性方法 等效线性方法是岩土地震工程中模拟波的传播的最常用的方法。 假定土体是粘弹性体，参照实验室得到的切线模量及阻尼比与剪应变幅值的关系曲线，对地震中每一单元的阻尼和模量重新赋值。Iterationtowardstrain-compatibleshearmodulusanddampingratio(afterKramer,1996)*ppt课件ppt课件*GeoHohai等效线性方法的特点 使用振动荷载的平均水平来估算每个单元的线性属性，并在振动过程中保持不变。在弱震阶段，单元会变得阻尼过大而刚度太小；在强震阶段，单元将会变得阻尼太小而刚度太大。对于不同部位不同运动水平的特性存在空间变异性。 不能计算永久变形。等效线性方法模型在加荷与卸荷时模量相同，不能计算土体在周期荷载作用下发生的剩余应变或位移。 塑形屈服模拟不合理。在塑性流动阶段，普遍认为应变增量张量是应力张量的函数，称之为“流动法则”。然而，等效线性方法使用的塑性理论认为应变张量（而不是应变增量张量）是应力张量的函数。因此，塑性屈服的模拟不合理。 大应变时误差大。等效线性方法所用割线模量在小应变时与非线性的切线模量很相近，但在大应变时二者相差很大，偏于不安全。 本构模型单一。等效线性方法本身的材料本构模型包括了应力应变的椭圆形方程，这种预设的方程形式减少了使用者的选择性，但却失去了选择其它形状的适用性。方法中使用迭代程序虽然部分考虑了不同的试验曲线形状，但是由于预先设定了模型形式，所以不能反映与频率无关的滞回圈。另外，模形是率无关的，因此不能考虑率相关性。*ppt课件ppt课件GeoHohai完全非线性分析方法FLAC3D采用完全非线性分析方法，基于显式差分方法，使用由周围区域真实密度得出的网格节点集中质量，求解全部运动方程。*ppt课件ppt课件GeoHohai完全非线性分析方法的特点 可以遵循任何指定的非线性本构模型。如果模型本身能够反映土体在动力作用下的滞回特性，则程序不需要另外提供阻尼参数。如果采用Rayleigh阻尼或局部（local）阻尼，则在动力计算中阻尼参数将保持不变。 采用非线性的材料定律，不同频率的波之间可以自然地出现干涉和混合，而等效线性方法做不到这一点。 由于采用了弹塑性模型，因此程序可以自动计算永久变形。 采用合理的塑性方程，使得塑性应变增量与应力相联系。 可以方便地进行不同本构模型的比较。 可以同时模拟压缩波和剪切波的传播及两者耦合作用时对材料的影响。在强震作用下，这种耦合作用的影响很重要，比如在摩擦型材料中，法向应力可能会动态地减小从而降低土体的抗剪强度。*ppt课件ppt课件GeoHohai使用弹塑性模型 附加考虑的因素： 阻尼，对于屈服面以下应力的循环 体积应变积累，是循环周数与幅值的函数 模量衰减，基于平均应变水平的最简单的弹塑性模型往往在描述累计塑性应变方面具有很好的效果，但是对于加速度放大系数的估算上效果不好。*ppt课件ppt课件GeoHohai弹塑性模型简单的理想弹塑性本构模型仅仅在发生屈服时才会出现滞回特性注意：即使这样粗糙的模型也能够作出连续的阻尼比和模量衰减曲线。在屈服条件下会产生体积改变，但通常都是剪胀。*ppt课件ppt课件*GeoHohaiF