(完整word版)用flac3d模拟基坑开挖
- 格式:doc
- 大小:73.50 KB
- 文档页数:12
运用flac3d模拟沟渠开挖前后岩体的应力应变报告一、问题的提出在天然山体与平原的的结合处开挖一条10m×4m×3m的沟渠,进行应力应变场的分析。
二、问题的分析与处理本模型的分析共分为三个阶段:第一阶段——前处理阶段:创建初始几何模型并划分网格,并在flac3d中显示网格体,然后定义材料的模型和材料参数,再定义边界条件和初始条件,这样便用flac3d建立起了计算的模型。
第二阶段——求解:在第一阶段建立好模型的基础上对模型进行变量的监控和求解运算。
第三阶段——后处理阶段:通过建立各种视口显示出各种运算的结果,以此判断模型的求解是否收敛。
对于沟渠的开挖模型及其运算,分三个阶段说明其生成及运算过程:1、前处理阶段由图可知,该模型可由两个radtunnel模型合并而成,分别以左侧模型p0(10,0,20)和右侧模型p0(10,10,10)为起始点。
步骤:1.打开flac3d软件,输入new命令全部清空原来数据,接着输入代码,如下:generate zone radtunnel p0 10 0 20 p1 0 0 20 p2 10 10 20 p3 10 0 0 p4 0 10 20 p5 10 10 0 p6 0 0 0 p7 0 10 0 p8 80 20 p9 10 0 7 p10 8 10 20 p11 10 10 7 p12 8 0 7 p13 8 10 7 generate zone radtunnel p0 10 10 10 p1 20 10 10 p2 10 0 10 p3 10 10 0 p4 20 0 10 p5 10 0 0 p6 20 10 0 p7 20 0 0 p8 12 10 10 p9 10 10 7 p10 12 0 10 p11 10 0 7 p12 12 10 7 p13 12 0 7plotadd surface blueshow生成网格模型如图1图1 几何网格模型2.定义材料的本构模型,本题为摩尔-库伦模型,然后对材料参数(体积模量,内摩擦角和剪胀角,粘聚力,抗拉强度)进行赋值,如没有输入其中一种参数,则系统默认为0。
1.FLAC3D知识基本介绍SimWe岩土工程结构的数值解是建立在满足基本方程(平衡方程、几何方程、本构方程)和边界条件下推导的。
由于基本方程和边界条件多以微分方程的形式出现,因此,将基本方程近假发改用差分方程(代数方程)表示,把求解微分方程的问题改换成求解代数方程的问题,这就是所谓的差分法。
差分法由来已久,但差分法需要求解高阶代数方程组,只有在计算机的出现,才使该法得以实施和发展。
FLAC3D(Fast Lagrangian Analysis of Continua)由美国Itasca公司开发的。
目前,FLAC 有二维和三维计算程序两个版本,二维计算程序V3.0以前的为DOS版本,V2.5版本仅仅能够使用计算机的基本内存(64K),所以,程序求解的最大结点数仅限于2000个以内。
1995年,FLAC2D已升级为V3.3的版本,其程序能够使用护展内存。
因此,大大发护展了计算规模。
FLAC3D是一个三维有限差分程序,目前已发展到V2.1版本。
FLAC3D的输入和一般的数值分析程序不同,它可以用交互的方式,从键盘输入各种命令,也可以写成命令(集)文件,类似于批处理,由文件来驱动。
因此,采用FLAC程序进行计算,必须了解各种命令关键词的功能,然后,按照计算顺序,将命令按先后,依次排列,形成可以完成一定计算任务的命令文件。
FLAC3D是二维的有限差分程序FLAC2D的护展,能够进行土质、岩石和其它材料的三维结构受力特性模拟和塑性流动分析。
调整三维网格中的多面体单元来拟合实际的结构。
单元材料可采用线性或非线性本构模型,在外力作用下,当材料发生屈服流动后,网格能够相应发变形和移动(大变形模式)。
FLAC3D采用的显式拉格朗日算法和混合-离散分区技术能够非常准确发模拟材料的塑性破坏和流动。
由于无须形成刚度矩阵,因此,基于较小内存空间就能够求解大范围的三维问题。
FLAC3D采用ANSI C++语言编写的。
FLAC3D有以下几个优点:1 对模拟塑性破坏和塑性流动采用的是“混合离散法“。
3D基坑开挖模拟学习目的通过2D网格扩展方式建立3D模型是GTS NX中比较特殊的创建3D模型的方法,不需要创建几何实体,仅需要划分2D网格后拉伸出3D网格,效率很高,特别适用于地层较为规则的模型。
通过本例题的学习,了解这种建模方法的要点和注意事项。
模型概况开挖深度:17m;开挖面积:72m×25m;支护形式:地连墙(30m)+2道砼撑+3道钢管撑;土层:共五层土,由上至下分别是杂填土、粘土、粉砂、中砂、强风化泥岩、中风化泥岩图1 模型示意图模型尺寸如下图所示:表2 结构参数表打开GTS NX主程序,点击左上角空白页新建模型文件。
确认为Z向,单位系统设定为KN、m、sec,初始参数按默认值,点击确定键。
3.导入几何线框导入准备好的基坑平面CAD图“基坑导入-3D.dwg”,选择导入dwg-3D,选择提供的CAD文件,点击确认或点击适用后关闭对话框(注意,只需点击一次,否则会导致重复导入)。
当然也可自行在GTS NX中绘制断面图,注意DWG文件要单独进行交叉分割,行成独立封闭面域,点击几何-交叉分割,选择所有线条,确认。
4.几何检查为了检查是否有重复的线段和未闭合的面域,采用检查重复、边-区域工具。
但是需要保证所有线段位于当前工作平面上(也就是格栅面),首先需要进行更改工作平面,确保XY为当前工作平面,左侧工作目录树中,双击XY即可,工作窗口中可以看到格栅面的变化。
端点以蓝色圆圈闪烁提示。
5.几何分组为了后期网格划分时便于快速找到相应几何对象,建议先进行分组,注意分组要在交叉分割之后完成,否则会打乱分组。
共分为四组,围檩、砼撑、钢管撑、土体。
其中围檩、土体比较简单,内撑的分组稍稍有些复杂,注意不要选错。
首先,左侧工作目录树中选中几何-右键-新几何组,创建3个新几个组,再选中该几何组,按F2(或是Fn+F2)进行重命名,注意改完名字后敲击回车。
然后,选中对应的几何组-右键-包含排除几何体,来选择需要分配的线条。
第43卷第6期 2009年6月上海交通大学学报J OU RNAL OF SHAN GHA I J IAO TON G UNIV ERSIT YVol.43No.6 J un.2009 收稿日期:2008207213基金项目:国家自然科学基金资助项目(50679041);上海市科学技术委员会资助项目(08201200903);江西省教育厅科学技术研究资助项目(G JJ 09367)作者简介:丁勇春(19792),男,江苏大丰市人,博士生,从事基坑与地下工程数值仿真方面的研究.王建华(联系人),男,教授,博士生导师,电话(Tel.):021*********,E 2mail :wjh417@. 文章编号:100622467(2009)0620976205基于FL AC3D 的基坑开挖与支护三维数值分析丁勇春1, 王建华1, 徐 斌1,2(1.上海交通大学土木工程系,上海200030;2.南昌工程学院土木工程系,南昌330029)摘 要:采用三维快速拉格朗日方法(FL AC3D )建立了考虑基坑分步开挖与支护全过程的三维动态计算模型,土体采用修正剑桥模型模拟,考虑了支护结构与土体的接触滑移作用,分析了基坑施工中围护墙变形、地表沉降、坑底隆起、坑外深层土体变形的基本特性.计算得到的地表沉降曲线与已有文献的经验沉降曲线基本一致,验证了计算结果的适用性.分析结果可为类似基坑工程的设计和施工提供有益参考.关键词:基坑开挖;数值分析;三维快速拉格朗日方法;修正剑桥模型;土与结构相互作用中图分类号:TU 473.2 文献标识码:AThree 2Dimensional Numerical Analysis of Braced ExcavationBased on F LAC 3DD I N G Yong 2chun 1, W A N G J i an 2hua 1, X U B i n1,2(1.Depart ment of Civil Engineering ,Shanghai Jiaotong University ,Shanghai 200030,China ;2.Depart ment of Civil Engineering ,Nanchang Instit ute of Technology ,Nanchang 330029,China )Abstract :Fast Lagrangian Analysis of Continua in 3Dimensions (FL AC3D )was employed to investigate t he deformatio n characteristics of a staged excavated and supported foundation pit.Modified Cam 2clay model was adopted to model t he soil behavior and t he interaction between soils and retaining st ruct ures was also taken into account in t he numerical model.The deflection of retaining walls ,t he settlement of gro und surface ,t he heave of excavation bottom ,and t he movement of t he deep st rata out side t he excavation were analyzed.The ground surface settlement curves of t he numerical model are basically consistent wit h t he empirical ones f rom t he existing literat ures ,so it s feasibility is proved.The numerical result s p rovide a usef ul reference for t he design and const ruction of similar deep excavation project s.Key words :braced excavation ;numerical analysis ;fast Lagrangian analysis of continua in 3dimensio ns (FL AC3D );modified Cam 2clay model ;soil 2st ruct ure interaction 基坑工程除了要保证基坑的整体性和支护结构的稳定安全外,还必须确定支护结构变形和基坑内外地层的变形,这样才能合理评估基坑施工对周围环境可能造成的不利影响,采取相应的工程措施,确保基坑施工的顺利进行,从而实现基于变形控制的基坑设计[1].目前基坑工程中常用的分析方法可以分为两类:行业和地方基坑设计规程所广为采用的竖向弹性地基梁法[2,3];基于连续介质力学的数值分析方法,如有限元法、有限差分法等.竖向弹性地基梁法只能对围护结构的内力和变形进行分析,不能计算支撑体系的内力和变形,更无法分析坑外地表沉降和坑底隆起变形,因此无法对基坑开挖引的环境影响进行分析和评价.数值分析方法可以弥补常规弹性地基梁法的不足.本文采用三维快速拉格朗日方法(FL AC3D )建立考虑支护结构与土体相互作用的基坑三维计算模型,土体采用修正剑桥模型,并考虑支护结构(围护墙和立柱桩)与土体间的接触问题,实现基坑分步开挖及支护的三维全过程动态分析.1 数值计算原理1.1 修正剑桥模型目前在岩土工程数值计算中常用的土体本构模型有两大类:一类是以Duncan 2Chang 模型为代表的非线性弹性模型;另一类是弹塑性模型,如Mohr 2Coulomb 模型、Drucker 2Prager 模型及修正剑桥模型等.Mohr 2Coulomb 模型和Drucker 2Prager 模型由于对加载和卸载采用同一模量,在得到合理围护结构侧向变形的同时,往往导致不合理的坑底隆起变形.修正剑桥模型能够反映土体加载与卸载模量的差异,考虑土体材料静水压力屈服特性和压硬性,在软土地基开挖分析中应用非常广泛.修正剑桥模型的屈服函数可表示为q 2p ′2+M 21-p ′0p ′=0(1)式中:q 为偏应力;p ′为平均有效应力;p ′0为先期固结压力;M 为p ′-q 平面内临界状态线(CSL )的坡度.修正剑桥模型屈服轨迹如图1所示.(a )p ′-q 平面内(b )主应力空间内图1 修正剑桥模型屈服轨迹Fig.1 Y ield locus of modified Cam 2clay model 修正剑桥模型的初始状态参数确定方法如下:(1)确定土体单元的竖向有效应力σ′v ;(2)根据现场静止侧压力系数K 0,确定土体单元侧向有效应力σ′h ,σ′h =K 0σ′v (2) (3)对每个土体单元,计算平均有效应力p ′和偏应力q :p ′=13(σ′1+σ′2+σ′3)q =12(σ′1-σ′2)2+(σ′2-σ′3)2+(σ′3-σ′1)2(3) (4)根据式(1)和超固结比OCR 计算先期固结压力p ′0,p ′0=p ′1+qM p ′2OCR (4) (5)确定正常固结线上p ′=1kPa 作用下单元的初始比体积υe ,0=υCSL +(λ-κ)ln 2(5)式中:λ为υ-ln p 平面内压缩线坡度;κ为υ-ln p 平面内回弹线坡度;υCSL 为临界状态线在p ′=1kPa 下的比体积.如根据土工试验能确定土样的不排水抗剪强度c u 和比体积υ0,则可根据下式确定参数υCSL[4],c u =M 2expυCSL -υ0λ(6) 如果已知土体颗粒比重G s 、土体液限w L 和土体塑性指数I p ,也可按下列经验公式确定参数[4]:υCSL =1+G s (w L +0.3I p )(7)1.2 土体与支护结构的接触采用数值方法进行基坑开挖与支护分析时涉及到是否需要考虑支护结构与土体间的相互接触问题.文献[5]中的研究结果表明,不考虑接触作用将导致计算结果与实际情况不符.主要原因可归结为:①无法反映支护结构与土体间由于接触面或土体的破坏引起的相互滑移;②结构与土体间存在拉应力,不能模拟结构与土体的脱离效应.FLAC3D 中的接触面单元具有单面特性,不同于G oodman 接触面单元.在每一计算时间步Δt 内,接触面节点与目标面之间的绝对法向侵入位移u n及相对剪切位移增量Δu s 均被计算,将其代入接触面本构方程就可以确定下一时刻(t +Δt )接触面上的法向力F (t +Δt )n和切向力F (t +Δt )s,F (t+Δt )n =k n u n A +σn A F (t+Δt )s=F (t )s+k s Δu (t+Δt/2)sA +τs A(8)式中:k n 为法向刚度;k s 为切向刚度;τs 为附加剪应力;A 为接触面节点的代表面积.接触单元服从库仑779 第6期丁勇春,等:基于FL AC3D 的基坑开挖与支护三维数值分析 剪切破坏准则和拉伸破坏准则.1.3 支护结构模型FLAC3D 提供了6种支护结构单元:梁单元、索单元、桩单元、壳单元、土工格栅单元和衬砌单元.对基坑而言,板桩式围护结构可采用衬砌单元模拟,桩列式和重力式挡土结构可采用实体模拟,水平内支撑如结构层楼板和混凝土或钢支撑可分别采用壳单元和梁单元模拟,坑外土钉和锚杆可采用索单元模拟,坑内立柱与立柱桩可采用桩单元模拟.2 计算模型及参数确定2.1 计算模型及边界条件分析模型取一方形基坑为研究对象,基坑平面尺寸为56m ×56m ,考虑模型的对称性后取1/4模型进行计算.基坑最大开挖深度H max =20m ,分5步开挖,每步开挖4m.1/4计算模型的三维尺寸为128m ×128m ×100m ,如图2所示.土体采用8节点6面体模拟,土体区段总数为13520,网格点总数为15309.模型外边界采用侧向约束,中心对称面采用对称边界,模型底部全约束.(a )整体模型(b )支护结构图2 计算模型Fig.2 Numerical model 计算中不考虑土体的分层以及基坑降水的影响,采用总应力法计算,相应的土体计算参数采用总应力指标.土体本构模型采用修正剑桥模型,不考虑地下水作用.土体重度γ=17.15kN/m 3,孔隙比e =1.2,侧压力系数K 0=0.5,κ=0.01,λ=0.14,M =1.2,泊松比μ=0.35.模型参数相当于上海软土地区第3层土体[6]的参数,按正常固结考虑(OCR =1).2.2 支护结构及接触面参数水平梁板支撑按刚度等效原则简化为壳单元,围护墙采用衬砌单元模拟,立柱和立柱桩采用桩单元模拟.结构单元总数为5175,结构单元节点数为2658.支护结构强度按C30混凝土考虑,考虑80%强度折减[7]后混凝土的弹性模量E c =24GPa ,μc =0.2,ρc =2500kg/m 3.坑边x =28m 和y =28m 处为两道连续墙,连续墙深度40m.坑内共设9根桩,桩直径0.8m ,间距8m ,桩长60m.共设5道水平支撑,支撑板厚0.12m.不同结构单元间采用共用节点法实现内力的传递.接触计算为物理非线性问题,有限元法在计算中往往存在收敛困难问题.FL AC3D 由于基于显式的全过程动力计算,因此,对不稳定问题不存在计算上的障碍.FLAC3D 中的衬砌单元与土体间的切向相互作用也具有单面特性,因而不能同时考虑围护墙与内外两侧土体的相互接触算法.在计算模型中采用以下近似处理办法:衬砌单元建立在墙外土体区域的外表面上以模拟围护墙与墙外土体的相互接触作用,围护墙与坑内土体的相互接触采用在坑内外土体间建立接触面单元,墙底处衬砌单元节点与坑内外土体网格点自由度耦合,假定墙底处结构单元节点与网格点变形协调.桩单元可直接实现桩与土界面的接触算法,通过加入桩端屈服弹簧可考虑桩的端承效应.参考文献[5],接触界面摩擦系数取为0.25,最大剪应力取为20k Pa.3 计算结果分析3.1 围护墙变形不同开挖深度(H )和不同平面位置(y )围护墙的侧向变形(δh )如图3所示.开挖至坑底后围护墙整体上抬约15.0mm ,墙体最大侧移-45.3mm,基坑开挖期间围护墙变形主要为侧向变形.(a )不同开挖深度 (b )不同平面位置图3 围护墙变形Fig.3 Deflection of retaining walls 由于第1道支撑的设置先于第1次4m 深土体的开挖,故第1次挖土时围护墙顶受到第1道水平支撑强大刚度的约束,墙顶位移几乎为零.后续阶段,墙体最大侧移随着开挖深度的增加而增大,同时最大侧移点位置亦随着开挖深度的增加而逐渐下879上 海 交 通 大 学 学 报第43卷 移,开挖至坑底后,墙身最大侧移点位于坑底开挖面附近.围护墙的最大侧移随着基坑开挖深度的增加呈非线性增长,但增长的幅度却有下降的趋势.围护墙的侧向变形具有明显的三维空间效应,基坑中心处(y =0)对称面墙体侧移最大,而基坑拐角位置(y =28m )墙体侧移最小,且为向坑内整体刚性平移.3.2 地表沉降不同开挖深度(H )和不同平面位置(y )坑外地表沉降(δv )如图4所示.图中,d 为计算点至坑边的距离.可见,地表未出现隆起,这是由于考虑了围护墙与土体间的接触特性,墙体与土体之间出现了滑移,墙体的上抬未引起地表的隆起.第1次挖土(H =4m )最大地表沉降点距坑边较近;随着挖土深度的不断增加,最大沉降点位置逐渐远离基坑,基坑开挖深度达12m 后,最大地表沉降点位置几乎不再变化.图4 地表沉降Fig.4 Ground surface settlement 地表沉降也具有空间效应,靠近基坑侧边中点(y =0)的地表沉降较大,最大地表沉降值为-36.3mm ;而靠近基坑拐角(y =28m )的地表沉降较小,最小地表沉降值仅为-16.4mm.地表沉降与围护墙的变形是相关联的,两者在空间上的变化规律保持一致. 图4(c )所示为中心对称面上不同开挖阶段相对地表沉降(δv /δv ,max )和相对距离(d/H )关系曲线与文献[8,9]中提出的经验沉降包络曲线的对比.可见,数值模拟结果与文献[9]中的更为吻合;d/H >2时本文与文献[8]中的相差较大.原因如下:文献[8]中假定地表沉降影响范围仅为2倍基坑开挖深度,而根据本文数值模拟结果和上海软土地区相关基坑的现场实测[10],上海软土地区地表沉降的影响范围一般为3~4倍基坑开挖深度.因此文献[8]中关于地表沉降影响范围的假定对上海软土地区的基坑而言显然偏小.由于计算模型的外截断边界采用侧向约束,边界面竖向可自由变形,引起了d/H >1.5范围计算所得相对地表沉降稍大于文献[9]中的包络曲线,这部分差异是由于计算模型的边界条件简化引起的.3.3 坑底隆起图5所示为不同开挖深度中心对称面(y =0)上坑底土体的隆起变形(δvb ).坑底土体的隆起主要由竖向开挖卸载效应引起,另外,墙体向基坑内的侧向变形也会进一步挤推坑内土体,造成坑底土体的回弹[11].坑底土体的隆起值随着开挖深度的增加呈非线性增长,但回弹增量有减小的趋势.图5 坑底隆起Fig.5 Heave of excavation bottom 开挖至坑底后坑底土体的竖向整体隆起量约160mm ,而围护墙的整体回弹仅15mm 左右,表明坑底土体与围护墙间发生了较大的相对滑移.如不考虑接触滑移作用,坑底回弹将引起过大的围护墙上抬,这显然不合理.另外,地表与围护墙相交处由于也考虑了接触作用,墙体的上抬才未引起相应位置地表的隆起,从而得到了与经验曲线基本吻合的地表沉降曲线.因此,在数值计算模型中有必要考虑支护结构与土体间的相互接触作用.3.4 坑外深层土体变形图6所示为开挖至坑底(H =20m )后中心对称面(y =0)上坑外深层土体竖向和侧向的位移场分布.可见,坑外浅层土体产生沉降,但最大沉降点位于距坑边一定距离处;坑外深层土体产生隆起,最大979 第6期丁勇春,等:基于FL AC3D 的基坑开挖与支护三维数值分析 隆起点位置位于土体与围护墙的接触处.邻近基坑围护墙处土体的侧向变形与墙体的变形相似,均为深层凸出型,最大侧向变形点深度位于基坑最终开挖面附近;随着至坑边距离的增加,土体最大侧向变形点位置逐步向地表面过渡,当距离超过一定范围,土体的最大侧移点位于地表面.(a )竖向变形(b )侧向变形图6 坑外深层土体变形Fig.6 Movement of strata outside excavation 对比图6中坑外土体的竖向变形与侧向变形可见,当基坑首层支撑刚度较大且先于浅层土体开挖而架设时,坑外浅层土体的变形以沉降为主,侧向变形较小.因此,对基坑周围建(构)筑物上部结构、浅埋地下市政管线及城市道路等而言,地表沉降和差异沉降是引起其损坏的主要原因;对基坑最终开挖面深度的坑外土体而言,侧向变形占主导,竖向变形较小,且存在隆起与下沉两种不同的形态.因此对深埋地下结构如高层建筑桩基础、城市高架桥梁桩基础、地铁区间隧道等而言,坑外深层土体的侧向变形是引起其损坏的主要原因.4 结 论(1)计算得到的坑外地表沉降曲线与已有文献的经验沉降包络曲线基本一致,验证了FL AC3D 在基坑开挖分析中的适用性和有效性.(2)基坑围护结构变形、地表沉降、坑底土体位移、坑外深层土体位移是相互关联的有机整体,基坑变形表现出明显的空间特性.(3)只要合理确定土体模型参数及土体与支护结构间的接触关系,三维数值分析在基坑工程设计施工方案优化和环境影响评估的应用中具有明显优越性.本文计算模型进行了一定简化,未考虑土体的分层,也未考虑降水引起的渗流和固结对基坑变形的影响,有待今后进一步研究.参考文献:[1] 刘建航,侯学渊.基坑工程手册[M ].北京:中国建筑工业出版社,1997.[2] J G J 120299,建筑基坑支护技术规程[S].[3] DBJ 08261297,基坑工程设计规程[S].[4] Wood D M.Soil behaviour and critical state soil me 2chanics [M ].Cambridge :Cambridge UniversityPress ,1990.[5] 范 巍,王建华,陈锦剑.连续墙与土体接触特性对深基坑变形分析的影响[J ].上海交通大学学报,2006,40(12):211822121.FAN Wei ,WAN G Jian 2hua ,CH EN Jin 2jian.The evaluation of deformation induced by excavation con 2sidering the properties of diaphragm 2soil interface [J ].Journal of Shangh ai Jiaotong U niversity ,2006,40(12):211822121.[6] D G J 0823722002,岩土工程勘察规范[S].[7] 谢百钩.粘土层开挖引致地盘移动之预测[D ].中国台湾:国立台湾科技大学,1999.[8] Clough G W ,O ′Rourke T D.Construction inducedmovements of in situ walls [C]//Proceedings ,ASCE Conference on Design and Perform ance of E arth R etai 2ning Structures .New Y ork :ASCE ,1990:4392470.[9] Hsieh P G ,Ou C Y.Shape of ground surface settle 2ment profiles caused by excavation [J ].C anadian G eotechnical Journal ,1998,35(6):100421017.[10] 丁勇春,王建华,徐中华,等.上海软土地区某深基坑施工监测分析[J ].西安建筑科技大学学报(自然科学版),2007,39(3):3332338.DIN G Y ong 2chun ,WAN G Jian 2hua ,XU Zhong 2hua ,et al .Monitoring analysis of a deep excavation inShanghai soft soil deposits [J ].Journal of Xi ’an U ni 2versity of Architecture &T echnology (N atural Science Edition),2007,39(3):3332338.[11] 刘国彬,黄院雄,侯学渊.基坑回弹的实用计算法[J ].土木工程学报,2000,33(4):61267.L IU Guo 2bin ,HUAN G Yuan 2xiong ,HOU Xue 2yuan.A practical method for calculating a heave of excavated foundation [J ].China Civil E ngineering Journal ,2000,33(4):61267.89上 海 交 通 大 学 学 报第43卷 。
FLAC3D简述与使用步骤FLAC3D是一种三维数值建模和数据分析软件,主要用于模拟和分析地下结构中的岩石和土壤行为。
它基于有限元方法,可以模拟地下开挖、地下水流、地震响应等复杂的地下工程问题,帮助工程师和地质学家做出准确的预测和决策。
在本文中,我们将对FLAC3D的概念和使用步骤进行简要介绍。
首先,我们来了解FLAC 3D的基本概念。
FLAC是Fast Lagrangian Analysis of Continua(快速拉格朗日连续体分析)的缩写,是一种用于建模和分析连续体力学问题的软件。
它采用了非线性弹性、塑性和损伤模型,并使用有限元离散化技术将复杂的问题转化为简单的网格模型。
FLAC 3D可以模拟岩土体的变形、破裂和失稳行为,帮助用户评估地下工程的安全性和可行性。
使用FLAC3D进行建模和分析的步骤如下:1.建立模型:在FLAC3D中,用户需要创建一个模型来描述地下结构。
模型可以包括岩石和土壤的几何形状、材料属性和边界条件等信息。
用户可以使用软件提供的几何建模工具创建模型,也可以导入其他CAD软件中的模型。
2.定义材料属性:在FLAC3D中,用户可以定义不同材料的物理和力学特性。
这些特性可以包括杨氏模量、泊松比、体积权重等。
用户可以根据实际材料的性质来设置这些参数,以便更真实地模拟地下结构的行为。
3.设置边界条件:在建模过程中,用户需要为模型设置适当的边界条件。
边界条件可以包括施加的加载、支撑结构和地下水流等。
用户可以通过定义加载的类型、大小和方向来模拟各种工程场景。
4.设定数值参数:在FLAC3D中,用户需要设置一些数值参数来控制数值计算的准确性和稳定性。
这些参数包括网格密度、时间步长和收敛准则等。
用户可以通过对不同参数的测试和调整来优化模拟结果的精度。
5.进行模拟和分析:完成模型设置后,用户可以运行FLAC3D来进行模拟和分析。
软件会根据用户定义的模型和参数对地下结构的行为进行预测和计算。
3.饱和土体中的开挖 3.1 问题陈述在不透水的基础饱和土体层中设计开挖。
土壤层1.2m 厚,潜水面水平是一个常数,与土壤顶部表面相一致。
开挖的洞穴将有一个8×8的正方形截面。
深5m 。
为本工程做准备,开挖地址已经被1m 厚的垂直封闭墙,该封闭墙深入到开挖洞穴下方2m 。
开挖以后,安装水泵用来降低开挖洞穴下部水平面。
问题是估计开挖洞穴底部由于开挖和排水引起的总体底鼓。
问题是根据对称性,分析时在三维空间取1/4部分考虑。
在土壤表面所在平面,开挖体中心为原点,定义直角坐标系统,z 轴方向指向下方。
在这个例子中,边界x =12m 和y =12m 作为对称面。
图3-1为本题的草图。
该土壤作为弹性材料。
土壤和水具有如下属性:弹性模量, K 390 MPa 剪切模量, G 280 MPa 土壤干密度, ρd 1200 kg/m 3水密度, ρw 1000 kg/m 3 墙密度, ρwall 1500 kg/m 3 渗透系数, k 10□12 m 2/Pa-s泊松比, n 0.3 流体体积模量, K f 2.0 GPa重力加速度g 约为10m/s 2。
在平衡状态的初始状态,各向同性应力为:'''4.0zzyy xx σσσ==3.2 建模过程FLAC 3D 模型尺寸为12m ×12m ×12m ,网格是尺寸为1m ×1m ×1m 的12×12×12的总体模型。
A fluid null mode l is assigned to the zones within the excavation and wall volumes , and a mechanical null model is assigned to the zones within the excavation. (See Figure 3.2 for a plot of the grid with the excavation removed.) 注意的是,在这个问题中,土壤饱和密度w d n ρρ+值与墙密度(wall density )相等。
;FLAC3D3.0在某隧道工程开挖支护中的应用;隧道建模命令流入下:newset log onset logfile yang.loggen zon radcyl p0 0 0 0 p1 9.0 0 0 p2 0 50 0 p3 0 0 8 &size 4 20 6 4 dim 6 5 6 5 rat 1 1 1 1 group 围岩gen zon cshell p0 0 0 0 p1 6.0 0 0 p2 0 50 0 p3 0 0 5.0 &size 4 20 6 4 dim 5.6 4.6 5.6 4.6 rat 1 1 1 1 group 初期支护gen zon cshell p0 0 0 0 p1 5.6 0 0 p2 0 50 0 p3 0 0 4.6 &size 4 20 6 4 dim 5.0 4.0 5.0 4.0 rat 1 1 1 1 group 二次衬砌 fill group 原岩gen zon radcyl p0 0 0 0 p1 0 0 -8.0 p2 0 50 0 p3 9.0 0 0 &size 4 20 6 4 dim 3 6 3 6 rat 1 1 1 1 group 围岩2gen zon cshell p0 0 0 0 p1 0 0 -3.0 p2 0 50 0 p3 6.0 0 0 &size 4 20 6 4 dim 2.6 5.6 2.6 5.6 rat 1 1 1 1 group 仰拱初期支护gen zon cshell p0 0 0 0 p1 0 0 -2.6 p2 0 50 0 p3 5.6 0 0 &size 4 20 6 4 dim 2 5 2 5 rat 1 1 1 1 group 仰拱二次衬砌 fill group 仰拱原岩gen zone reflect normal -1 0 0gen zone radtun p0 0 0 0 p1 45 0 0 p2 0 50 0 p3 0 0 20 &size 3 20 3 12 dim 9 8 9 8 rat 1 1 1 1.1 group 围岩3gen zon reflect dip 0 ori 0 0 0 range x 0 9 y 0 50 z 8 20gen zon reflect dip 0 ori 0 0 0 range x 9 45 y 0 50 z 0 20gen zon reflect dip 90 dd 270 ori 0 0 0 range x 0 9 y 0 50 z 8 20gen zon reflect dip 90 dd 270 ori 0 0 0 range x 0 9 y 0 50 z -8 -20gen zon reflect dip 90 dd 270 ori 0 0 0 range x 9 45 y 0 50 z -20 20gen zon brick p0 -45 0 -20 p1 -45 0 -40 p2 -45 50 -20 p3 45 0 -20 &size 5 20 6 rat 1.1 1 1 group 围岩4save tun_model.sav;假设围岩岩体符合mohr-coulomb本构模型,给围岩赋参数命令流如下,; mohr-coulomb modelmodel mohrdef derives_mod1=E_mod1/(2.0*(1.0+p_ratio1))b_mod1=E_mod1/(3.0*(1.0-2.0*p_ratio1))s_mod2=E_mod2/(2.0*(1.0+p_ratio2))b_mod2=E_mod2/(3.0*(1.0-2.0*p_ratio2))endset E_mod1=0.6e9 p_ratio1=0.27 E_mod2=0.8e9 p_ratio2=0.26deriveprop bulk b_mod1 shear s_mod1 cohe 1.8e6 tens 0.8e6 fric 30 range z 4.5 20 prop bulk b_mod2 shear s_mod2 cohe 2.8e6 tens 1.0e6 fric 35 range z -40 4.5 ini dens=2300set grav 0 0 -10; boundary and initial conditionsapply szz -1.4e6 range z 19.9 20.1fix z range z -40.1 -39.1fix x range x -45.1 -44.9fix x range x 44.9 45.1fix y range y 49.9 50.1hist unbalhist gp xdis 6.0,0,0hist gp zdis 0,0,5hist gp xdis 6.0,50,0hist gp zdis 0,50,5plot hist 3solvesave tun_nature.sav;对后面计算而言,模型建立时岩体在开挖前认为位移已经终了,因此需要对位移进行“清零”,而应力可以保留。
new;网格建立;;;;;;;;;;;;;;;;;;;;;;;;;gen zone brick p0 90 0 -30 p1 202 0 -30 p2 90 4 -30 p3 90 0 0 size 112 4 30 ratio 1 1 1gen zone brick p0 90 0 -30 p1 90 0 0 p2 90 4 -30 p3 0 0 -30 size 30 4 25 ratio 1 1 1.1gen zone brick p0 90 0 -30 p1 0 0 -30 p2 90 4 -30 p3 90 0 -75 size 25 4 18 ratio 1.1 1 1.1gen zone brick p0 90 0 -30 p1 90 0 -75 p2 90 4 -30 p3 202 0 -30 size 18 4 112 ratio 1.1 1 1 gen zone brick p0 202 0 -30 p1 292 0 -30 p2 202 4 -30 p3 202 0 0 size 25 4 30 ratio 1.1 1 1 gen zone brick p0 202 0 -30 p1 202 0 -75 p2 202 4 -30 p3 292 0 -30 size 18 4 25 ratio 1.1 1 1.1;分组;;;;;;;;;;;;;;;;;;;;;;;;;;group 1 range x 90 110 y 0 4 z -30 0group 1 range x 180 202 y 0 4 z -30 0group 2 range group 1 not;建立连续墙单元;;;;;;;;;;;;;;;;;;;;;;;;;;gen separate 1gen merge 1e-4 range x 90 110 y 0 4 z -30.1 -29.9gen merge 1e-4 range x 180 202 y 0 4 z -30.1 -29.9attach face range x 89.99 90.01 y 0.0 4.0 z -29.9 0attach face range x 109.99 110.01 y 0.0 4.0 z -29.9 0attach face range x 179.99 180.01 y 0.0 4.0 z -29.9 0attach face range x 201.99 202.01 y 0.0 4.0 z -29.9 0sel liner id 1 crossdiag group 2 range x 89.9 90.1 y -0.1 4.1 z -30.1 0.1sel liner id 2 crossdiag group 2 range x 109.9 110.1 y -0.1 4.1 z -30.1 0.1sel liner id 3 crossdiag group 2 range x 179.9 180.1 y -0.1 4.1 z -30.1 0.1sel liner id 4 crossdiag group 2 range x 201.9 202.1 y -0.1 4.1 z -30.1 0.1sel liner id 1 prop isotropic (2.0e10, 0.20) thickness 0.8 density 2.5e3 &cs_nk=4e9 cs_sk=4e9 &cs_ncut=4e7 cs_scoh=4e7 cs_scohres=0 cs_sfric=20.0 &range x 89.9 90.1 y -0.1 4.1 z -30.1 0.1sel liner id 2 prop isotropic (2.0e10, 0.20) thickness 0.8 density 2.5e3 &cs_nk=4e9 cs_sk=4e9 &cs_ncut=4e7 cs_scoh=4e7 cs_scohres=0 cs_sfric=20.0 &range x 109.9 110.1 y -0.1 4.1 z -30.1 0.1sel liner id 3 prop isotropic (2.0e10, 0.20) thickness 0.8 density 2.5e3 &cs_nk=4e9 cs_sk=4e9 &cs_ncut=4e7 cs_scoh=4e7 cs_scohres=0 cs_sfric=20.0 &range x 179.9 180.1 y -0.1 4.1 z -30.1 0.1sel liner id 4 prop isotropic (2.0e10, 0.20) thickness 0.8 density 2.5e3 &cs_nk=4e9 cs_sk=4e9 &cs_ncut=4e7 cs_scoh=4e7 cs_scohres=0 cs_sfric=20.0 &range x 201.9 202.1 y -0.1 4.1 z -30.1 0.1;定义支撑结构;;;;;;;;;;;;;;;;;;;;;;;def struct_install1loop i(1,3)structx_zz=-1.0*5.0*(i-1)structx_xx0=90.0structx_xx1=110.0structx_yy=2.0commandsel beam id=2 begin (structx_xx0,structx_yy,structx_zz) end (structx_xx1,structx_yy,structx_zz) nseg=10sel beam id=2 prop dens=0.000 emod=1.0e-6 nu=0.2 &xcarea=0.80 xcj=10.94e-2 xciy=6.67e-2 xciz=4.27e-2 ydirection=(0 0 -1) ;1000x800endcommandendloopendstruct_install1def struct_install2loop i(1,3)structx_zz=-1.0*5.0*(i-1)structx_xx0=180.0structx_xx1=202.0structx_yy=2.0commandsel beam id=3 begin (structx_xx0,structx_yy,structx_zz) end (structx_xx1,structx_yy,structx_zz) nseg=11sel beam id=3 prop dens=0.000 emod=1.0e-6 nu=0.2 &xcarea=0.80 xcj=10.94e-2 xciy=6.67e-2 xciz=4.27e-2 ydirection=(0 0 -1) ;1000x800endcommandendloopendstruct_install2;建立结构单元分组;;;;;;;;;;;;;;;;;;;;;;;;;;;sel group linerwall range sel linersel group struct1 range sel beam x (90.0 110.0) z (-0.1 0.1)sel group struct2 range sel beam x (90.0 110.0) z (-5.1 -4.9)sel group struct3 range sel beam x (90.0 110.0) z (-10.1 -9.9)sel group struct4 range sel beam x (180.0 202.0) z (-0.1 0.1)sel group struct5 range sel beam x (180.0 202.0) z (-5.1 -4.9)sel group struct6 range sel beam x (180.0 202.0) z (-10.1 -9.9);删除beam单元的linksel dele link range sel beam z (-30 0);建立liner间的节点间的刚性linkdef merge_link0node_num=0node_pnt0 = nd_headloop while node_pnt0 # null ;寻找总节点数,注:不能自己任生成node,程序缺省的方式为连续生成无不连续node_num = node_num+1node_pnt0 = nd_next(node_pnt0)endloopnode_num_minus1 = node_num-1link_id=30000loop ii (1,node_num_minus1)node_pnt1 = nd_find(ii)xxa = nd_pos(node_pnt1,2,1)yya = nd_pos(node_pnt1,2,2)zza = nd_pos(node_pnt1,2,3)ii_plus1 = ii+1loop jj (ii_plus1,node_num)node_pnt2 = nd_find(jj)xxb = nd_pos(node_pnt2,2,1)yyb = nd_pos(node_pnt2,2,2)zzb = nd_pos(node_pnt2,2,3)node_dist = sqrt((xxa-xxb)^2+(yya-yyb)^2+(zza-zzb)^2)dist_tol = 1e-1if node_dist <= dist_tol thenlink_pnt1 = nd_link(node_pnt1)link_pnt2 = nd_link(node_pnt2);if link_pnt1 # null then; temp1 = lk_delete(link_pnt1);endifif link_pnt2 # null thentemp2 = lk_delete(link_pnt2)endiflink_id = link_id+1command ;生成新link(6自由度全固结),大的node的id作为target node,小的node的id作为source node,需注意不同情况下的灵活调整sel set link node_tol=dist_tolsel link id=link_id jj target = node tgt_num =ii ;指定link的ID;sel link ii target = node tgt_num = jj ;不指定link的id,自动生成sel link attach xdir=rigid ydir=rigid zdir=rigid xrdir=rigid yrdir=rigid zrdir=rigid range id=link_idendcommandendifendloopendloopendmerge_link0;设置土层材料参数;;;;;;;;;;;;;;;;;;;;;;;;;;;;def b_s_modb_mod =e_mod/(3.0*(1.0-2.0*p_ratio))s_mod =e_mod/(2.0*(1.0+p_ratio))endmodel elasticset e_mod 100e6set p_ratio 0.3b_s_modprop bu=b_mod sh=s_modini dens 1800 range z -75 0def ini_szzszz0=0szzgrad=1800*10commandini szz add szz0 grad 0 0 szzgrad range z -75 0endcommandendini_szzdef ini_sxx_syypnt=zone_headloop while pnt # nullval=k0*z_szz(pnt)z_sxx(pnt)=valz_syy(pnt)=valpnt=z_next(pnt)endloopendset k0=0.50ini_sxx_syy;定义边界处的结构边界条件;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;cyc 0sel node local xdir=(0,1,0) ydir=(0,0,1) range x 89.9 90.1 y -0.1 4.1 z -30.1 0.1 sel node local xdir=(0,1,0) ydir=(0,0,-1) range x 109.9 110.1 y -0.1 4.1 z -30.1 0.1 sel node local xdir=(0,1,0) ydir=(0,0,1) range x 179.9 180.1 y -0.1 4.1 z -30.1 0.1 sel node local xdir=(0,1,0) ydir=(0,0,-1) range x 201.9 202.1 y -0.1 4.1 z -30.1 0.1sel node fix lsys range x 89.9 90.1 y -0.1 0.1 z -30.1 0.1sel node fix lsys range x 89.9 90.1 y 3.9 4.1 z -30.1 0.1sel node fix lsys range x 109.9 110.1 y -0.1 0.1 z -30.1 0.1sel node fix lsys range x 109.9 110.1 y 3.9 4.1 z -30.1 0.1sel node fix lsys range x 179.9 180.1 y -0.1 0.1 z -30.1 0.1sel node fix lsys range x 179.9 180.1 y 3.9 4.1 z -30.1 0.1sel node fix lsys range x 201.9 202.1 y -0.1 0.1 z -30.1 0.1sel node fix lsys range x 201.9 202.1 y 3.9 4.1 z -30.1 0.1sel node fix x yr zr range x 89.9 90.1 y -0.1 0.1 z -30.1 0.1sel node fix x yr zr range x 89.9 90.1 y 3.9 4.1 z -30.1 0.1sel node fix x yr zr range x 109.9 110.1 y -0.1 0.1 z -30.1 0.1sel node fix x yr zr range x 109.9 110.1 y 3.9 4.1 z -30.1 0.1sel node fix x yr zr range x 179.9 180.1 y -0.1 0.1 z -30.1 0.1sel node fix x yr zr range x 179.9 180.1 y 3.9 4.1 z -30.1 0.1sel node fix x yr zr range x 201.9 202.1 y -0.1 0.1 z -30.1 0.1sel node fix x yr zr range x 201.9 202.1 y 3.9 4.1 z -30.1 0.1sel node fix y range x 89.9 90.1 y 0.0 4.0 z -0.1 0.1sel node fix y range x 109.9 110.1 y 0.0 4.0 z -0.1 0.1sel node fix y range x 179.9 180.1 y 0.0 4.0 z -0.1 0.1sel node fix y range x 201.9 202.1 y 0.0 4.0 z -0.1 0.1;set plot meta;plot set rot 20 0 30 ba wh color=on cent=(10 20 0) mag=3.81;set outp node_local_sys.wmf;plot add sel geom black red link=off node=off id=off shrink=0 scale=0.03 nodesys=on range group linerwall any group struct1 any;pl ha;固定边界条件;;;;;;;;;;;;;;;;;;;;;;;;;;fix x range x -0.1 0.1fix x range x 291.9 292.1fix y range y -0.1 0.1fix y range y 3.9 4.1fix x y z range z -75.1 -74.9set grav 0,0,-10solvesave elas.sav;删除侧面内外土体间的连接约束;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;attach delete range x 89.99 90.01 y 0.0 4.0 z -29.9 0attach delete range x 109.99 110.01 y 0.0 4.0 z -29.9 0attach delete range x 179.99 180.01 y 0.0 4.0 z -29.9 0attach delete range x 201.99 202.01 y 0.0 4.0 z -29.9 0;在墙内土体的外侧建立接触面;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;interface 1 face range group 1 x 89.99 90.01 y 0.0 4.0 z -29.9 0interface 2 face range group 1 x 109.99 110.01 y 0.0 4.0 z -29.9 0interface 3 face range group 1 x 179.99 180.01 y 0.0 4.0 z -29.9 0interface 4 face range group 1 x 201.99 202.01 y 0.0 4.0 z -29.9 0interface 1 prop kn=4e8 ks=4e8 tens=5e3 coh=0.0 fric=20 ;接触面参数interface 2 prop kn=4e8 ks=4e8 tens=5e3 coh=0.0 fric=20 ;接触面参数interface 3 prop kn=4e8 ks=4e8 tens=5e3 coh=0.0 fric=20 ;接触面参数interface 4 prop kn=4e8 ks=4e8 tens=5e3 coh=0.0 fric=20 ;接触面参数interface 1 maxedge=1interface 2 maxedge=1interface 3 maxedge=1interface 4 maxedge=1;interface 1 prop kn=4e8 ks=4e8 tens=1e10 sbratio=100;plot set ba wh;pl ske interface red blue attach cyan green;set outp interface_attachment.wmf;pl ha;重新定义连续墙参数;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;sel liner id 1 prop isotropic (2.0e10, 0.20) &cs_nk=4e9 cs_sk=4e9 &cs_scoh=4e7 cs_scohres=0.0 cs_sfric=0.0 range x 89.9 90.1 y -0.1 4.1 z -30.1 0.1sel liner id 2 prop isotropic (2.0e10, 0.20) &cs_nk=4e9 cs_sk=4e9 &cs_scoh=4e7 cs_scohres=0.0 cs_sfric=0.0 range x 109.9 110.1 y-0.1 4.1 z -30.1 0.1sel liner id 3 prop isotropic (2.0e10, 0.20) &cs_nk=4e9 cs_sk=4e9 &cs_scoh=4e7 cs_scohres=0.0 cs_sfric=0.0 range x 179.9 180.1 y -0.1 4.1 z -30.1 0.1sel liner id 4 prop isotropic (2.0e10, 0.20) &cs_nk=4e9 cs_sk=4e9 &cs_scoh=4e7 cs_scohres=0.0 cs_sfric=0.0 range x 201.9 202.1 y -0.1 4.1 z -30.1 0.1;重新定义墙底约束条件;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;def redef_wall_end_link1node_pnt = nd_headlink_id=100000loop while node_pnt # nullnode_id = nd_id(node_pnt)xx = nd_pos(node_pnt,2,1)yy = nd_pos(node_pnt,2,2)zz = nd_pos(node_pnt,2,3)link_pnt = nd_link(node_pnt)dist_x = sqrt((xx-90.0)^2+(zz+30.0)^2)if dist_x <=dist_tol thenif link_pnt # null thentemp1 = lk_delete(link_pnt)\link_id = link_id+1commandsel set link node_tol = dist_tolsel link id=link_id node_id target zonesel link attach xdir=rigid ydir=rigid zdir=rigid xrdir=free yrdir=free zrdir=free range id=link_idendcommandendifendifnode_pnt = nd_next(node_pnt)endloopendredef_wall_end_link1def redef_wall_end_link2node_pnt = nd_headlink_id=150000loop while node_pnt # nullnode_id = nd_id(node_pnt)xx = nd_pos(node_pnt,2,1)yy = nd_pos(node_pnt,2,2)zz = nd_pos(node_pnt,2,3)link_pnt = nd_link(node_pnt)dist_x = sqrt((xx-110.0)^2+(zz+30.0)^2)dist_tol = 1e-1if dist_x <=dist_tol thenif link_pnt # null thenif yy < 85.0temp1 = lk_delete(link_pnt)\link_id = link_id+1commandsel set link node_tol = dist_tolsel link id=link_id node_id target zonesel link attach xdir=rigid ydir=rigid zdir=rigid xrdir=free yrdir=free zrdir=free range id=link_idendcommandendifendifendifnode_pnt = nd_next(node_pnt)endloopendredef_wall_end_link2def redef_wall_end_link3node_pnt = nd_headlink_id=200000loop while node_pnt # nullnode_id = nd_id(node_pnt)xx = nd_pos(node_pnt,2,1)yy = nd_pos(node_pnt,2,2)zz = nd_pos(node_pnt,2,3)link_pnt = nd_link(node_pnt)dist_x = sqrt((xx-180.0)^2+(zz+30.0)^2)dist_tol = 1e-1if dist_x <=dist_tol thenif link_pnt # null thentemp1 = lk_delete(link_pnt)\link_id = link_id+1commandsel set link node_tol = dist_tolsel link id=link_id node_id target zonesel link attach xdir=rigid ydir=rigid zdir=rigid xrdir=free yrdir=free zrdir=free range id=link_idendcommandendifendifnode_pnt = nd_next(node_pnt)endloopendredef_wall_end_link3def redef_wall_end_link4node_pnt = nd_headlink_id=250000loop while node_pnt # nullnode_id = nd_id(node_pnt)xx = nd_pos(node_pnt,2,1)yy = nd_pos(node_pnt,2,2)zz = nd_pos(node_pnt,2,3)link_pnt = nd_link(node_pnt)dist_x = sqrt((xx-202.0)^2+(zz+30.0)^2)dist_tol = 1e-1if dist_x <=dist_tol thenif link_pnt # null thentemp1 = lk_delete(link_pnt)\link_id = link_id+1commandsel set link node_tol = dist_tolsel link id=link_id node_id target zonesel link attach xdir=rigid ydir=rigid zdir=rigid xrdir=free yrdir=free zrdir=free range id=link_idendcommandendifendifnode_pnt = nd_next(node_pnt)endloopendredef_wall_end_link4;剑桥模型;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;model cam-clay;cam-clay模型则不需定义弹性模量(E、G、K)等参数,自动计算;cam-clay模型中需确定8个模型参数(①-⑧),手册property中的初始比体积cv(v0)和shear 无须给定def install_proppnt=zone_headloop while pnt # nullabs_sxx=abs(z_sxx(pnt)) ;|sxx|abs_syy=abs(z_syy(pnt)) ;|syy|abs_szz=abs(z_szz(pnt)) ;|szz|p0=(abs_sxx+abs_syy+abs_szz)/3.0;cam-clay模型中p、q均须为正值,p0由初应力场确定,故cam-clam定义模型参数前须先已知初应力p0_effective=p0-z_pp(pnt) ;p0';q0=sqrt(((abs_sxx-abs_syy)^2+(abs_syy-abs_szz)^2+(abs_szz-abs_sxx)^2)*0.5)q0=sqrt(((abs_sxx-abs_syy)^2+(abs_syy-abs_szz)^2+(abs_szz-abs_sxx)^2)*0.5+3.0*((z_sxy(pnt)) ^2+(z_sxz(pnt))^2+(z_syz(pnt))^2))z_prop(pnt,'mm')=6.0*sin(fai*degrad)/(3.0-sin(fai*degrad)) ;①注三角函数中需将角度转化为弧度temp1=q0/(z_prop(pnt,'mm')*p0_effective)pc0=p0_effective*(1.0+temp1^2)*OCR ;先期有效固结压力,用于确定屈服面v0=1.0+_e0z_prop(pnt,'cam_cp')=p0_effective ;★重要参数,否则不能正确计算有效应力,提示出错"Mean effective pressure is negative"z_prop(pnt,'mpc')=pc0 ;②z_prop(pnt,'poisson')=p_ratio ;③z_prop(pnt,'lambda')=_lambda ;④z_prop(pnt,'kappa')=_kappa ;⑤z_prop(pnt,'mp1')=_mp1 ;⑥z_prop(pnt,'mv_l')=v0+_lambda*ln(2.0*_cu/(z_prop(pnt,'mm')*_mp1))+(_lambda-_kappa)*l n(2.0) ;⑦z_prop(pnt,'bulk_bound')=100*40e6 ;⑧;z_prop(pnt,'bulk_bound')=100*(s_mod+4.0/3.0*s_mod) ;弹性体模上界Kmax;自动确定Kmax时会出现“property bad”错误提示;因为弹性上界对计算结果无影响,在不提示Kmax太小的性况下,取值越小计算收敛越快pnt=z_next(pnt)endloopendset p_ratio=0.25 fai=34.5 _lambda=0.14 _kappa=0.012 _mp1=1e3 _e0=1.2 _cu=10e3 OCR=1.0 ;模型所需参数install_propsolvesave model_cam.sav。