FLAC-FLAC3D基础与应用-(2)课件

FLAC/FLAC3D基础与应用Email:2011年11月24日1几个问题什么是FLAC？为什么要用FLAC？FLAC能做什么？FLAC为何这么流行？怎么学FLAC？2什么是FLAC？FastLagrangianAnalysisofContinua3为什么要用FLAC？Abaqus： 745条Plaxis： 80条ADINA： 310条4FLAC能做什么？岩土工程中的绝大多数问题土力学、岩石力学、防灾减灾、隧道、地下空间等采矿工程中的大部分问题水工结构中的部分问题结构工程貌似无能为力国际通用的岩土工程专业分析程序5FLAC为何这么流行？CharlesFairhurst美国工程院、瑞典皇家工程院院士，国际岩石力学学科和岩石力学学会创始人之一，历任国际岩石力学学会主席和副主席，国际岩石力学学会Muller奖、美国岩石力学学会终生成就奖获得者。PeterCundall美国工程院、英国皇家工程院院士，国际资深计算岩石力学学家。6怎么学习FLAC？问：好学不？答：好！问：难不？答：入门容易深入较难7报告安排第一讲：基本介绍、静力分析、前后处理第二讲：接触面、FISH语言、流固耦合分析第三讲：动力分析、自定义本构、结构单元第四讲：FLAC（2D）基本介绍与应用实例讨论8第一讲FLAC3D基本介绍、静力分析、前后处理9FLAC大应变、小应变计算模式.丰富的本构模型、提供自定义的本构模型功能接触面可以模拟不同材料的接触流固耦合实现土体的固结与渗流拥有各种功能的结构单元类型，模拟土与结构的相互作用强大的动力分析功能.流变分析，拥有粘弹性模型和粘塑性模型热力学分析.Shearstrainratecontours10FLAC3D

与FLAC类似，是FLAC的三维版本与FLAC拥有相同的优点upstreamdownstream

11可选模块可选模块包括:热力学,热-力学耦合,热-流体-力学耦合包括热传导和对流;粘弹,粘朔性(蠕变)材料模型;动力学分析,并可以模拟静边界和自由域使用C++定义自己的模型核废料储存中的热力学研究问题12NewFeaturesinFLAC3DVersion4.0模拟颗粒状材料的硬化模型自动网格重画功能，解决bad-geometry问题.改进的interface更快的渗流计算更新的动力计算功能13Lagrangian格式动量平衡方程F(t)m牛顿运动定律对于连续体在静力平衡条件下，加速度项为0，方程变为平衡方程14自由落体的模拟G=mgS=1/2gt2=20m命令流：configdyngenzonbrisize111inixmul0.1ym0.1zm0.1modelelaspropbulk3e8shear1e8inidens1000setgrav00-10solveage215自由落体的模拟(movie)16FLAC3D的求解过程平衡方程(动量方程)应力—应变关系(本构模型)Gauss定律单元积分应变率速度节点力新的应力对所有的网格节点对所有单元17FLAC3D中的本构模型开挖模型null3个弹性模型各向同性弹性横观各向同性弹性正交各向同性弹性8个弹塑性模型Drucker-Prager模型、Morh-Coulomb模型、应变硬化/软化模型、遍布节理模型、双线性应变硬化/软化遍布节理模型、修正剑桥模型和胡克布朗模型18分析问题的过程建立网格初始条件边界条件初始应力平衡外荷载求解前处理后处理19初始应力的生成为什么要单独列出？分析过程中出现的很多问题都与初始应力是否合理有关手册中的例子五花八门是所有后续分析的基础！！生成方法弹性求解更改强度参数的弹塑性求解设置初始应力的弹塑性求解存在水压力的初始应力生成水下建筑的初始应力生成20弹性求解模型尺寸单元数量密度KGu1×1×2(m3)1×1×2200030MPa10MPa0.35genzonbrisize112melaspropbulk3e7shear1e7fixzranz0fixxranx0fixxranx1fixyrany0fixyrany1inidens2000setgrav00-10solveStep=162sz=-40e3sx=-21.54e321更改强度参数的弹塑性求解模型尺寸单元数量密度KGcfu1×1×2(m3)1×1×2200030MPa10MPa10kPa150.35genzonbrisize112modelmohrpropbulk3e7shear1e7c1e10f15ten1e10fixzranz0fixxranx0fixxranx1fixyrany0fixyrany1inidens2000setgrav00-10solvepropbulk3e7shear1e7c10e3f15ten0solveStep=163sz=-40e3sx=-21.54e3Or:solveelastic22设置初始应力的弹塑性求解模型尺寸单元数量密度KGcfu1×1×2(m3)1×1×2200030MPa10MPa10kPa150.35genzonbrisize112modelmohrpropbulk3e7shear1e7c10e3f15ten0fixzranz0fixxranx0fixxranx1fixyrany0fixyrany1inidens2000iniszz-40e3grad0020e3ranz02inisyy-20e3grad0010e3ranz02inisxx-20e3grad0010e3ranz02setgrav00-10solveStep=0sz=-40e3sx=-20e323存在水压力的初始应力生成(1)模型尺寸单元数量饱和密度KGcfu水位线孔隙率1×1×2(m3)1×1×2200030MPa10MPa10kPa150.351m0.5genzonebricksize112modelmohrpropbulk3e7shear1e7coh10e3fri15ten0fixzranz0fixxranx0fixxranx1fixyrany0fixyrany1inidens2000ranz01inidens1500ranz12iniszz-35e3grad0020e3ranz01inisyy-22.5e3grad0015e3ranz01inisxx-22.5e3grad0015e3ranz01iniszz-30e3grad0015e3ranz12inisyy-15e3grad007.5e3ranz12inisxx-15e3grad007.5e3ranz12inipp10e3grad00-10e3ranz01setgrav00-10solveStep=0sz=-35e3sx=-22.5e3rd=rs–n×s×rf24存在水压力的初始应力生成(2)模型尺寸单元数量饱和密度KGcfu水位线孔隙率1×1×2(m3)1×1×2200030MPa10MPa10kPa150.351m0.5configfluidgenzonbrisize112modelelaspropbu3e7sh1e7inidens1500modelfl_isoinifdens=1000fmod0proppor0.5setgrav00-10waterdens1500watertableface001,011,111,101inipp10e3grad00-10e3ranz01iniszz-30e3grad0015e3ranz12...setfluidoffsolveStep=142sz=-40e3sx=-24e3单元数较少产生的误差25水下建筑的初始应力生成模型尺寸单元数量饱和密度KGcfu水位线1×1×2(m3)1×1×2200030MPa10MPa10kPa150.353mgenzonbrisize112modelmpropbulk3e7shear1e7c10e10f15ten1e10fixzranz0fixxranx0fixxranx1fixyrany0fixyrany1inidens2000ranz02iniszz-50e3grad0020e3ranz02inisyy-40e3grad0015e3ranz02inisxx-40e3grad0015e3ranz02inipp30e3grad00-10e3ranz02appnstress-10e3ranz2setgrav00-10solveStep=0sz=-50e3sx=-40e326前后处理27FLAC3D的前处理FLAC3D网格生成的关键特征：FLAC3D是命令驱动.使用FLAC3D内置基元进行形状组合可形成复杂网格.用户自定义FISH函数可以用来修改基元网格以创建更为复杂的网格.第三方软件导入.28FLAC3D网格基元块体退化块体楔体金字塔四面体柱体径向块体径向隧道径向柱体柱状壳体柱状交叉隧道交叉29Tecplot的后处理30Tecplot的云图和等值线Tecplot中云图的效果与flac3d的shadeon的效果差不多,但是其出三维等值线的功能是卓越的.31Tecplot的切片功能Tecplot的slice切片功能.和CAD一样,可以任意切剖面出图.最大优点是可以几个剖面同时出图.32整体和截面的数据以单元形式存在列表中,可自由选择和组合多个单元出图.33外界模型的导入复杂网格的生成难度大接口编写不同软件之间的精度差异Group的定义采用.flac3d的文件格式.flac3d文件的格式G11.0e+001.0e+001.0e+00ZB812345678ZGROUPSoil123impgrid&expgrid仅限于网格支持abaqus,ansys,HM,etc.3wzones3sec34复杂网格的检查FLAC3D生成的复杂网格attachface无接触面时检查整体模型是否存在sub-grid有接触面时给定范围进行检查其他软件导入的复杂模型网格划分的检查genmerge弹性模型modelelastic求解“独立”节点“畸形”单元35模型的检查FLAC本身的Check功能十分有限错误提示很少十分开放的工作平台(自由落体)检查的基本步骤网格检查(如前所述)边界条件检查速度约束条件 plogpfixredsk荷载条件 plofapredsk模型检查模型赋值 ploblockmodel参数赋值 ploblockprop***初始应力检查(如后所述)36初始应力检查办法初始应力的计算时间不会“非常长”经常检查模型的响应ploconszz(syy,sxx) 应力场ploconzdis(ydis,xdis) 位移场ploblosta 屈服状态plogpfixredsk 速度约束条件plofapredsk 体力plohist(unbal) 不平衡力plointerfacenstress(sstress) 接触面单元37第二讲FLAC3D接触面、FISH语言、流固耦合分析FISH语言FISH语言简介软件自带的编程语言是否一定要学？视情况而定，需要时查询FISH变量即可语法简单xxx……end_xxx注意事项与FLAC本身的关键字冲突保留字不可缩写变量可不定义，因此注意检查程序printfishtable,extra等命令使用一个最简单的FISH程序defabcabc=1+2*3abcd=1.0/2.0endabcprintfish数据格式函数与变量都可以在FISH函数中进行赋值，赋值操作与常规的编程语言类似，按照运算符的优先级先后顺序来执行。函数和变量的赋值遵守数据类型的规则，即整型的计算结果为整型，浮点型的计算结果为浮点型，因此读者在进行除法运算、开方运算时都需要将数据类型设置为浮点型，数字尽量使用小数点以保证运算正确。变量和函数名的命名规则是不能以数字开头，不能含有中文，并且不能包含如下的字符。

.,*/+-ˆ=<>#()[]@;’"变量和函数名不能与FLAC3D、FISH的保留字相冲突，不要采用过于简单的单词，比如a，hist等，这些都与保留字相冲突。即使程序中存在与保留字相冲突的变量，FLAC3D也不会提供任何提示，所以提醒读者在编制FISH程序时尽量使用较长的、复杂的变量和函数名。函数与变量对变量进行赋值时，不能使用当前函数的函数名放在“=”的右边，比如采用下面的定义

abcd=abc+1.0

在FISH程序执行时会提出错误，因为这样会形成递归调用，这种调用方式在FISH程序中是不允许的。变量和函数的作用是全局的，在命令中的任何地方修改变量的值都会立即生效，因此在实际应用中尽量避免不同的函数中含有相同的变量，因为这样可能会造成赋值错误，并难以检查。在FLAC3D中可以用如下的命令来引用FISH函数和变量PRINT用于查看函数和变量的数值；HISTORY命令可以对函数和变量的数值进行记录；SET命令用于变量的赋值。主要语句选择语句CASEOF

表达式…默认语句CASEn1…表达式的值为n1时的语句CASEn2…表达式的值为n2时的语句ENDCASE主要语句条件语句IF

条件表达式[THEN]…[ELSE]…ENDIFFISH中条件运算符没有“并”、“或”、“否”这样的符号ifaa>1.0ifaa<2.0

执行语句

endifendif表达“1<aa<2”的条件主要语句循环语句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单元遍历语句节点遍历语句主要语句命令语句COMMAND…（FLAC3D命令）ENDCOMMAMD应用实例：土体的模量随小主应力变化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应用实例获得最大位移的大小及发生位置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_value

maxdisp_value=disp_gpmaxdisp_gpid=gp_id(p_gp)endifp_gp=gp_next(p_gp)endloopendfind_max_dispprintmaxdisp_valuemaxdisp_gpidFISH的编写习惯第一步第二步defabcendabcdefabc

p_gp=gp_headloopwhilep_gp#nullp_gp=gp_next(p_gp)endloopendabc第三步第四步defabcp_gp=gp_headloopwhilep_gp#null

commandendcommandp_gp=gp_next(p_gp)endloopendabcdefabcp_gp=gp_headloopwhilep_gp#nullcommand

appnstress

…

endcommandp_gp=gp_next(p_gp)endloopendabcFISH与建模两个圆形隧道的连接部分变直径的隧道部分FISH检查主要采用PRINTfish命令查看变量的赋值是否合理，主要检查值为0的函数和变量，因为FISH程序中一般定义的变量都有实际的意义，输出为0的变量很可能是与保留字相冲突的变量（如a就是apply的保留字）或者由于编写笔误产生的变量（如数字0与大写字母O）。流固耦合分析很难！流-固耦合分析(单相流)基本功能理论框架计算模式渗流边界条件，初始条件单渗流计算及渗流耦合计算基本功能渗流各向同性、各向异性不同的渗流模型和属性流体压力，涌入量，渗漏量和不渗水边界抽水井、点源、体积源饱和渗流可采用显式差分法、隐式差分法非饱和渗流采用显式差分法渗流-固体-热的耦合流体和固体的耦合程度依赖于土体颗粒(骨架)的压缩程度，用Biot系数表示颗粒的可压缩程度。循环荷载引起的动水压力变化和土体液化理论框架准静态Biot理论多孔介质中遵循Darcy定律的单相渗流描述多孔介质中流体渗流的变量孔隙水压力，饱和度，特定排水向量的三个分量质量守恒定律达西定律本构定律考虑流体响应孔隙水压力改变，饱和度改变，体积应变改变和温度改变有效应力计算无渗流模式不设置CONFIGFluid孔压不改变设置孔压分布INITIALppWATERtableWATERdensityWATERtablefaceSETgravity手动设置干湿密度渗流模式设置CONFIGfluid设置土体干密度

渗流模型MODELfl_isotropicMODELfl_anisotropicMODELfl_nullSETfluidoffsetWATERbulk=0渗流边界条件，初始条件默认的边界条件是不透水边界孔隙压力自由(不透水边界)固定孔隙水压力(透水边界)如：水井

孔隙压力，孔隙率，饱和度和流体属性的初始分布可以用INITIAL命令或者PROPERTY命令定义。单渗流计算及渗流耦合计算时间比例完全耦合分析方法孔压固定分析(有效应力分析)单渗流得到孔压分布无渗流计算——孔压的力学响应流-固耦合计算单渗流得到孔压分布用途：排水沟；抽水井；耦合计算计算步骤CONFIGfluidSETmechoffSETfluidimpliciton/offMODELfl_;PROPSTEP;SOLVEage;SETfluidratioSETfluidoffmechonPROPbiot_c0(orINIfmod0)无渗流计算——孔压的力学响应不排水短期响应两种分析方法：干法和湿法干法：Ku=K+a2M两种破坏形式WATER或INI获得常孔压，不排水的c,φ(孔压改变较小)φ=0,c=cu(M>>K+4/3G)湿法：耦合体系的短期行为使用排水的K,c,φ若SETfluidoff,Biot_mod(fmod)真实力学过程的特征时间流体扩散过程的特征时间完全耦合分析方法时间比例短期行为(不排水)ts(分析时间)<<tc(耦合扩散时间)忽略渗流影响长期行为(排水)ts>>tc施加扰动的属性流体扰动：渗流可不与力学过程耦合力学扰动：耦合等级取决于流固刚度比流固刚度比流-固耦合计算CONFIGfluid;M(Kf);K(渗透系数)真实，则FLAC3D默认耦合计算Δp→Δev→sΔev→Δp预估流/力特征时间耦合计算前先达到一个平衡状态SETfluidonmechoff;SETfluidoffmechon;STEPSETmechforce;SETmechsubstepnauto;SETfluidsubstepm(=1)STEP:渗流步足够小流固耦合的计算方法手动调整的STEP求解

SETfluidonmechoff STEP SETfluidoffmechon STEP 主从进程的SOLVE求解SETmechforceSETmechsubstepnauto (从进程)SETmechsubstepm (主进程)SOLVEage自动STEP求解STEP渗流问题(CONFIGfluid)分析步骤时间比例(ts,tc)稳态不排水状态相当扰动类型力学扰动孔压扰动流固刚度比Rk是否>>>1完全耦合模式时间比例相当；力学扰动心墙土坝的渗流(1)newconfigfluidsetfluidoffgenzonbrickp000-10size20110genzonbrickp0500p11500p2510p3905p41510p5915p61105p71115size1015groupsoilgroupdamranx57z-50groupdamranid201aid211aid221aid231aid241agroupdamranid202aid212aid222aid232aid242amepropbu3e7sh1e7inipp0grad00-10e3ranz0-10inidens2000modelfl_isoproppor0.5perm1e-10inifden1000ften-1e10inisat0.0ranz05modelfl_nullrangrodam;inipp0rangrodamfixzranz-10fixxranx0fixxranx20fixysetgrav10solvesaveelastic.sav网格模型初始孔压心墙土坝的渗流(2)restelastic.savinixd0yd0zd0xv0yv0zv0appnstress-40e3grad0010e3ranz04x09solvesavepressure.sav竖向应力沉降心墙土坝的渗流(3)restpressure.savsetfluidonmechoffinifmod2e3ften0.0rangrosoilinixd0yd0zd0xv0yv0zv0apppp40e3grad00-10e3ranz04x09apppp0ranz0x1520histid=10zoneppid215solveRatio=1Ratio=1E-5荷载引起的地基土体的超孔隙水压力

10m3mK=500MPa,G=300MPa,c=10kPaf=1540kPa20m计算文件;---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步的求解。计算结果对主从进程法的讨论SETmechforce

设置一个不平衡力的大小，达到这个不平衡力系统认为暂时达到平衡状态；SETmechsubstepnauto设置力学进程为从进程，在主进程每执行一步中必须执行n步，当系统达到平衡时也可以少于n步；SETfluidsubstepm设置流体进程为主进程。讨论1：对收敛准则进行对比分析

不平衡力（force）：1E3、5E3、1E4、5E4不平衡力比（ratio）：1E-4、1E-3、1E-2、1E-1采用1E-3的收敛准则既可以满足计算流固耦合过程中的计算精度要求，同时又具有较高的计算效率对主从进程法的讨论(2)讨论2：子步数的影响

setmechsub100fluidsub10setmechsub10fluidsub10setmechsub1fluidsub1设置合理的子步数也很重要。设置过大，则会导致计算时间大大增加，过小又会造成计算结果的误差。真空预压的简单模拟孔压边界条件ts>>tc长期分析(排水)Rk>>1骨架很软孔压扰动进行biot_mod调整砂层软土层粘土层PVD2m8m10mDatafile:数值分析过程(movie)FLAC3D非线性动力分析非常复杂！SaidbyProf.PeterCundall76为什么要用FLAC做动力分析？FLAC可以模拟体系（土，岩石，结构，流体）受到的外部动力荷载（比如地震）或内部动力荷载（比如基础振动、爆炸）。可以计算塑性引起的永久变形以及孔隙水压力的消散。土动力学中常用的等效线性方法无法直接处理上述问题。77动力模拟的3个重要问题动力荷载与边界条件材料响应与阻尼土体液化78动力荷载动力输入的类型加速度时程速度时程应力(压力)时程力时程APPLYINTERIOR(内部)TABLEFISH79Quiet边界静态(quiet,粘性)边界LysmerandKuhlemeyer(1969)模型边界法向和切向设置独立的阻尼器性能对于法向p波和s波能很好的吸收对于倾斜入射的波和Rayleigh波也有所吸收，但存在反射人工边界仍应当足够远80Quiet边界应用内部振动(如隧道中的列车振动问题)

动力荷载直接施加在节点上使用Quiet边界减小人工边界上的反射不需要FF边界外部荷载的底部边界

软土地基上的地震荷载不适合用加速度或速度边界条件使用应力条件t=-2Csrvs地震底部输入的侧向边界

扭曲了入射波quietquietquiet81Free-field边界Cundalletal.(1980)自由场网格与主体网格的耦合粘性阻尼器，自由场网格的不平衡力施加到主体网格边界上设置条件底部水平，重力方向为z向侧面垂直，法向分别为x,y向其他边界条件在APPLYff之前相当于一个阻尼器82Free-field边界APPLYff将边界上单元的属性、条件和变量全部转移ff单元上；设置以后主体网格上的改动将不会被FF边界所响应可存在任意的本构模型以及流体耦合(仅竖向)FF边界进行小变形计算，主体网格可大变形，FF边界上的变形要相对较小存在attach的边界将不能设置FF边界边界上的Interface将不能连续动力边界设置需在FF边界设置之前83Free-field边界与动力荷载模型底部边界fix——施加速度或加速度荷载——刚性边界Free——施加应力时程荷载——柔性边界对于软弱的地基不适合施加速度（加速度荷载），而应当施加应力荷载Notethatthereisafactorof2becausetheinputenergydividesintoadownward-&upward-propagatingwave.842.材料响应与阻尼连续的非线性，表观模量随着应变的增大而降低对所有循环应变等级均存在滞回特性，因此导致随着循环应变的增加阻尼比增大。阻尼是率相关的。对于复杂波形的各个成分都产生阻尼。剪切应变会产生的体积应变，相应的，随着剪应变循环次数的增加体积应变逐渐积累。85材料响应土体在循环荷载作用下呈现出模量衰减和能量消散的特点，那么如何用非线性数值方法对其进行模拟呢？Nonlinearcharacteristicsofsoils(MartinandSeed,1979)86试验得到的阻尼比、割线模量随循环剪应变的曲线

0.00010.0010.010.11Shear

Strain

Amplitude

(%)01020304050Damping

Ratio

(%)0.00010.0010.010.11Shear

Strain

Amplitude

(%)0.00.10.20.30.40.50.60.70.80.91.0Normalized

Shear

Modulus,

G/GmaxMid-RangeSandCurve(Seed&Idriss,1970)SandFillInland:Friction=32,hr=0.47,Go=440SandFillunderRockDike:Friction=30,hr=0.43,Go=44087等效线性方法等效线性方法是岩土地震工程中模拟波的传播的最常用的方法。假定土体是粘弹性体，参照实验室得到的切线模量及阻尼比与剪应变幅值的关系曲线，对地震中每一单元的阻尼和模量重新赋值。Iterationtowardstrain-compatibleshearmodulusanddampingratio(afterKramer,1996)88等效线性方法的特点使用振动荷载的平均水平来估算每个单元的线性属性，并在振动过程中保持不变。在弱震阶段，单元会变得阻尼过大而刚度太小；在强震阶段，单元将会变得阻尼太小而刚度太大。对于不同部位不同运动水平的特性存在空间变异性。不能计算永久变形。等效线性方法模型在加荷与卸荷时模量相同，不能计算土体在周期荷载作用下发生的剩余应变或位移。塑形屈服模拟不合理。在塑性流动阶段，普遍认为应变增量张量是应力张量的函数，称之为“流动法则”。然而，等效线性方法使用的塑性理论认为应变张量（而不是应变增量张量）是应力张量的函数。因此，塑性屈服的模拟不合理。大应变时误差大。等效线性方法所用割线模量在小应变时与非线性的切线模量很相近，但在大应变时二者相差很大，偏于不安全。本构模型单一。等效线性方法本身的材料本构模型包括了应力应变的椭圆形方程，这种预设的方程形式减少了使用者的选择性，但却失去了选择其它形状的适用性。方法中使用迭代程序虽然部分考虑了不同的试验曲线形状，但是由于预先设定了模型形式，所以不能反映与频率无关的滞回圈。另外，模形是率无关的，因此不能考虑率相关性。89完全非线性分析方法FLAC3D采用完全非线性分析方法，基于显式差分方法，使用由周围区域真实密度得出的网格节点集中质量，求解全部运动方程。90完全非线性分析方法的特点可以遵循任何指定的非线性本构模型。如果模型本身能够反映土体在动力作用下的滞回特性，则程序不需要另外提供阻尼参数。如果采用Rayleigh阻尼或局部（local）阻尼，则在动力计算中阻尼参数将保持不变。采用非线性的材料定律，不同频率的波之间可以自然地出现干涉和混合，而等效线性方法做不到这一点。由于采用了弹塑性模型，因此程序可以自动计算永久变形。采用合理的塑性方程，使得塑性应变增量与应力相联系。可以方便地进行不同本构模型的比较。可以同时模拟压缩波和剪切波的传播及两者耦合作用时对材料的影响。在强震作用下，这种耦合作用的影响很重要，比如在摩擦型材料中，法向应力可能会动态地减小从而降低土体的抗剪强度。91使用弹塑性模型附加考虑的因素：阻尼，对于屈服面以下应力的循环体积应变积累，是循环周数与幅值的函数模量衰减，基于平均应变水平的表格最简单的弹塑性模型往往在描述累计塑性应变方面具有很好的效果，但是对于加速度放大系数的估算上效果不好。92弹塑性模型简单的理想弹塑性本构模型仅仅在发生屈服时才会出现滞回特性strainstress注意：即使这样粗糙的模型也能够作出连续的阻尼比和模量衰减曲线。在屈服条件下会产生体积改变，但通常都是剪胀。93FLAC3D中的阻尼比弹塑性本构模型使用瑞利（粘性）阻尼弹塑性本构模型使用滞后（HD）阻尼.复杂本构模型拥有连续的屈服应力应变关系和对应的加卸载响应（Wang，UBCSand）94瑞利阻尼瑞利阻尼最初应用于结构和弹性体的动力计算中，以减弱系统的自然振动模式的振幅。在计算时，假设动力方程中的阻尼矩阵C与刚度矩阵K和质量矩阵M有关：瑞利阻尼中的质量分量相当于连接每个节点和地面的阻尼器，而刚度分量则相当于连接单元之间的阻尼器。虽然两个阻尼器本身是与频率有关的，但是通过选取合适的系数，可以在有限的频率范围内近似获得频率无关的响应。95瑞利阻尼frequencyratioofdampingtocritical混合仅有刚度分量仅有质量分量采用叠加的方法得到的阻尼比在较大的频率范围内保持定值（3:1）Combinedcurvereachesminimumat:96中心频率的选择

假设弹性模型计算土石坝不同材料的功率谱曲线(assumingelasticmaterial)frequencyfrequencyfrequencyfrequencyfrequency97阻尼比的选择

根据65%最大应变来选择阻尼比和模量衰减参参数98瑞利阻尼的缺点必须指定中心频率，有时须在一些充满矛盾的数据当中进行选择（场地响应或地震平均频率）刚度比例项会导致时间步减小99滞后阻尼FLACVersion5.0andFLAC3DVersion3.0提供了滞后阻尼功能，该项阻尼独立于材料模型之外，实质是包含了一个与应变幅值相关的切线剪切模量乘子如果割线模量是由衰减曲线来确定的，那么切线模量可以表示为：FromSeed&Idriss(1970)Go=small-strainshearmodulus100滞后阻尼表观应变是偏应变，根据先前反转的点进行累计得到的。使用堆栈的数据结构可以描述应变反转现象。因此，在“迷你滞回圈”中的能量损失与主滞回圈的计算是一样的。elasticmodelwithhysteresisdamping特定的曲线可以使用滞后阻尼来描述，并在本构模型计算中根据应变的大小获得不同的切线模量。101滞后阻尼的特点可以直接采用动力试验中的模量衰减曲线；相对于瑞利阻尼而言，滞后阻尼不影响动力计算的时间步；可以应用于任意的材料模型，且可以与其它阻尼格式同时使用。一个缺点：模量衰减曲线一致，但阻尼比曲线存在差别“Good”fittoSeed&IdrisdataforG/Gmax(sigmoidal3-parameterfunction)–noteinconsistentdampingresult.G/GmaxD-%ofcritical102滞后阻尼的比较滞后的默认模型滞后的三参数模型滞后的四参数模型滞后的哈丁模型滞后阻尼低循环应变下得到的阻尼比要小于试验结果，这会导致低级的噪声，尤其在高频情况下。可以在中心频率上增加一个小量的Rayleigh阻尼(~0.2%刚度比例)，这样也不会降低时步；若初始应力不为0，剪应力-剪应变曲线可能不匹配。因此在生成初始应力时就要调用Hyst阻尼；Hyst阻尼不仅会增加能量损失，还会导致在大循环应变下的平均剪切模量的降低，在输入波的基频接近共振频率的时候，由于可能会导致动力反应幅值的增大；Hyst阻尼之前要做一次弹性无阻尼求解，以获得发生循环应变的最大水平，若循环应变过大导致剪切模量过多的降低，那么用Hyst阻尼是有问题的；即使应变较小，使用屈服模型也会增大应变，因此若有广泛屈服的现象，则使用屈服模型，不用Hyst阻尼104动孔压的生成——液化干沙剪应变循环加载试验初始加载阶段，沙土通常先压实再膨胀。卸载时，沙土遵循与加载相似的路径，但在零应变时，有些残余体积应变存在。取决于初始孔隙率，这可能代表纯粹的压实假定孔隙中充满水对于常体积测试，有效应力降低，孔隙水压保持不变对于常荷载测试，(例如，盒子上法向荷载固定)，孔隙水压增加，有效应力减小有效应力为零时发生液化105动孔压的生成——液化因此孔隙水压增加不是液化的基本原因由于颗粒间(重组以后)的低接触力导致有效应力的减小描述液化的模型高级模型：BSHP(边界面低塑性本构模型,Wangetal.1990)简单模型：MC+体积应变增量模型Finn模型：Byrne模型：106UBCTOT:Post-liquefactionAssignpost-liquefactionpropertiestgt107PL-Finn开始初始液化判断?零有效应力计算YES非零有效应力计算零有效应力判断?NO普通的Finn模型YESNO动力时间？结束NOYES0：零有效应力状态1：非零有效应力状态0有效应力非0有效应力108地震波的调整基线校正对于地震分析的加速度时程，其积分得到的速度和位移应归0美国地质调查研究所BasicStrong-MotionAccelerogramProcessingSoftware(BAP)对网格施加一个固定速度从而使残余的位移变为0动力荷载的频率与单元尺寸的双向调整高频的输入要求单元尺寸很小一定的单元尺寸对应输入的最大频率一般进行滤波处理滤掉低能量的高频FFT.FISOriginSeismoSignal109地震波的调整El-Centro波FFT修正后的时程修正后FFT5Hz110PL-Finn模型的应用阪神地震的码头分析沉箱顶部水平残余位移最大达5m，平均为3.5m，残余沉降为1~2m，海侧倾斜角3~5o。

seaSiltreplacedsandstonelandSandcaissonclayFLAC3Dgrid8050zones10386gps可液化的砂土111分析参数GroupnameConstitutivemodelrd(kg/m3)E(MPa)mc(kPa)j(o)clayMC1350500.333020seasiltMC1250200.33030replacedsandMC1350150.33037landsandMC135013.70.33036stoneMC15501000.33040caissonElastic350020000.17‐‐GroupnameFluidmodelK

(cm/s)PorositynDampingratioDliquefiedparametersPL-Finnparametersclayfl_iso1.0E-60.450.05‐‐seasiltfl_iso1.0E-050.450.05‐‐replacedsandfl_iso1.0E-030.450.05ByrneModelDr=40%C1=0.751C2=0.533C3=

0pprc=0.99k0=3105.4n0=0.3225k1=5503.1n1=0.1739landsandfl_iso1.0E-030.450.05ByrneModelDr=25%C1=2.432C2=0.164C3=

0pprc=0.99k0=3105.4n0=0.3225k1=5503.1n1=0.1739stonefl_iso1.0E-010.450.05‐‐caissonfl_null‐‐0.05‐‐基本力学参数流体参数动力参数PL-Finn模型参数112输入加速度南北向地震波最大值=0.6g竖直向地震波最大值=0.2g113接触面与自由边界条件模型周围自由场边界条件114计算结果位移震后水平位移等值线沉箱顶点处的水平位移、沉降时程曲线3.44m1.99m115计算结果（续）超孔压比震后的超孔压比云图超孔压比时程曲线A:置换砂B:海砂（底部）C:海砂（中部）116小结FLAC的动力反应分析非常复杂，用户在进行动力分析之前要具有熟练的静力分析基础FLAC的动力计算需要花费很多时间，主要原因是由于有限差分法本身需要较小的时间步，因此建议在进行动力分析前首先进行简单模型的分析“小例子”虽然会额外花费一定时间，但对于熟悉边界条件、荷载施加、阻尼选择等内容来说非常必要自定义本构模型118自定义本构模型的基本方法必要性试验总结的本构模型特定条件下的本构模型交叉学科的本构模型二次开发环境自定义本构模型的功能自定义本构模型的基本方法119二次开发环境FLAC3D采用面向对象的语言标准C++编写本构模型都是以动态连接库文件(.DLL文件)的形式提供VC++2005或更高版本的开发环境优点自定义的本构模型和软件自带的本构模型的执行效率处在同一个水平自定义本构模型(.DLL文件)适用于高版本的FLAC(2D)、3DEC、UDEC等其他Itasca软件中120自定义本构模型的功能主要功能：对给出的应变增量得到新的应力辅助功能：模型名称、版本读写操作模型文件的编写基类(classConstitutiveModel)的描述成员函数的描述模型的注册模型与FLAC3D之间的信息交换模型状态指示器的描述121自定义本构模型的基本方法头文件(usermodel.h)中进行新的本构模型派生类的声明修改模型的ID(>100)、名称和版本修改派生类的私有成员C++文件(usermodel.cpp)中修改模型结构(UserModel::UserModel(boolbRegister):ConstitutiveModel)constchar**UserModel::Properties()函数模型的参数名称字符串constchar**UserModel::States()函数计算过程中的状态指示器122自定义本构模型的基本方法doubleUserModel::GetProperty()和voidUserModel::SetProperty()函数constchar*UserModel::Initialize()函数参数和状态指示器的初始化，并对派生类声明中定义的私有变量进行赋值constchar*UserModel::Run()

函数由应变增量计算得到应力增量，从而获得新的应力constchar*UserModel::SaveRestore()函数对计算结果进行保存。程序的调试在VC++的工程设置中将FLAC3D软件中的EXE文件路径加入到程序的调试范围中，并将FLAC3D自带的DLL文件加入到附加动态链接库(AdditionalDLLs)中，然后在Initialize()或Run()函数中设置断点，进行调试；在程序文件中加入return()语句，这样可以将希望得到的变量值以错误提示的形式在FLAC3D窗口中得到。123相关文件mohr.hmohr.cpp124一个例子(Duncan-Chang)125其他成功的例子南京水科院双屈服面模型双屈服面中心截面沉降云图水平位移云图小主应力云图126FLAC基本介绍与应用实例127FLACGIIC基本操作Interface的建模流体计算128GIIC基本操作差分网格的建模思路成层土建模养成Mark的好习惯修剪模型撤销的最好方法模型草图计算网格数量确定模型生成方法不规则形状生成Attach,interface参数化运行FLAC129Build选项卡130Alter131Material132In-situ133StructureBeam，梁单元 ——可以用于各种类型的支护模拟，包括开挖支护、隧道中的支架等。可以在其两侧连接Interface单元来模拟岩土介质中的挡土墙，还可以通过Interface单元与FLAC网格相连以模拟土工格栅。Liner，衬砌单元 ——主要用于隧道衬砌，包括混凝土初衬或喷射混凝土初衬。Cable，锚索单元 ——不能承受弯矩，可以施加预应力，常用于模拟受拉构件，包括岩石中的锚杆。Pile，