蠕变本构
(Viscous)
Model visc（H-N）
一. 力-位移定律
Maxwell体－M体
F
u uH uN
总位移 弹性部 分 粘性部 分
uH
uN
将虎克体和牛顿体的本构关系带入做变换可得Visc本构 力-位移定律：
FБайду номын сангаас
H-N体示意图
Model Visc（H-N）
二. 中心差分方案（将应力张量做球应力张量和偏应力张量分解-增量 关系）
返回目录
自定义本构模型的功能
主要功能：对给出的应变增量得到新的应力 辅助功能： 模型名称、版本 读写操作 模型文件的编写 基类(class Constitutive Model)的描述 成员函数的描述 模型的注册 模型与FLAC3D之间的信息交换 模型状态指示器的描述
bool bFast=ps>stnS.Resoltopris(&dPrinMin,&dPrinMid,&dPrinMax,&a Dir,uDim, &icase, &sdif, &psdif);
double dPrinMinCopy = dPrinMin; double dPrinMidCopy = dPrinMid; double dPrinMaxCopy = dPrinMax; /* --- Mohr-Coulomb failure criterion --- */ double dFsurf = dPrinMin - dNPH * dPrinMax + dCSN; /* --- Tensile failure criteria --- */ double dTsurf = dTension - dPrinMax; double dPdiv = -dTsurf + (dPrinMin - dNPH * dTension + dCSN) * dBISC;
/* --- tests for failure */ if (dFsurf < 0.0 && dPdiv < 0.0) { iPlas = 1; /* --- shear failure: correction to ps->incips->l stresses ---*/ ps->mState = (unsigned long)(ps->mState | 0x01); dPrinMin -= dFsurf * dSC1; dPrinMid -= dFsurf * dSC2; dPrinMax -= dFsurf * dSC3; } else if (dTsurf < 0.0 && dPdiv > 0.0) { iPlas = 2; /* --- tension failure: correction to ps->incips->l stresses ---*/ ps->mState = (unsigned long)(ps->mState | 0x02); double dTco = dE21 * dTsurf; dPrinMin += dTco; dPrinMid += dTco; dPrinMax = dTension; } if (iPlas) { ps->stnS.Resoltoglob(dPrinMin,dPrinMid, dPrinMax, aDir, dPrinMinCopy,dPrinMidCopy,dPrinMaxCopy, uDim, icase, sdif, psdif, bFast); ps->bViscous = false; // Inhibit stiffness-damping terms } else { ps->bViscous = true; // Allow stiffness-damping terms } return(0); }

常用模型信息传递指针变量
返回目录
模型源程序分解
返回目录
静力本构
(Mohr-Coulomb)
MC本构
① ②
①－剪切屈服 ②－拉伸屈服
1. 屈服函数
剪切屈服
2. 塑性势函数
非关联
拉屈服
关联
剪切屈服修正
写成线性函数S：
拉伸屈服修正
写成线性函数S：
const char *UserMohrModel::Initialize(unsigned uDim,State *) { if ((uDim!=2)&&(uDim!=3)) return("Illegal dimension in UserMohr constitutive model"); dE1 = dBulk + d4d3 * dShear; dE2 = dBulk - d2d3 * dShear; dG2 = 2.0 * dShear; double dRsin = sin(dFriction * dDegRad); dNPH = (1.0 + dRsin) / (1.0 - dRsin); dCSN = 2.0 * dCohesion * sqrt(dNPH); if (dFriction) { double dApex = dCohesion * cos(dFriction * dDegRad) / dRsin; dTension = dTension < dApex ? dTension : dApex; } dRsin = sin(dDilation * dDegRad); dRnps = (1.0 + dRsin) / (1.0 - dRsin); double dRa = dE1 - dRnps * dE2; double dRb = dE2 - dRnps * dE1; double dRd = dRa - dRb * dNPH; dSC1 = dRa / dRd; dSC3 = dRb / dRd; dSC2 = dE2 * (1.0 - dRnps) / dRd; dBISC = sqrt(1.0 + dNPH * dNPH) + dNPH; dE21 = dE2 / dE1; return(0); }
ITASCA(武汉)咨询有限公司
报告人: 朱永生 2007.04.28 依泰斯卡(武汉)
内

容
必要性
试验总结的本构模型 特定条件下的本构模型 交叉学科的本构模型


二次开发环境
自定义本构模型的功能 自定义本构模型的基本方法 常用模型信息传递指针变量 模型源程序分解
二次开发环境
其中:
d d 偏应力张量 和偏应变张量 的本构 关 式中：
总应力：
球应力张量 关系
iso
和球应变张量 kk 的本构
Model visc（H-N）
三. 源程序分解
const char *UserViscousModel::Initialize(unsigned,State *) { G 求解系数：dGD2V= ； dG2 = 2.0 * dShear; 2 if (dViscosity <= 0.0) dGD2V = 0.0; 1. 如果粘滞系数 <=0， dGD2V=0； else dGD2V = 0.5 * dShear / dViscosity; 2. 否则dGD2V为真值。 return(0); } const char *UserViscousModel::Run(unsigned uDim,State *ps) { if ((uDim!=3)&&(uDim!=2)) return("Illegal dimension in UserViscousModel"); double dD = dGD2V * (ps->bCreep ? ps->dTimeStep : 0.0); if (dD > 0.5) return("Timestep too large for UserViscousModel"); 求解系数：dD= 1. 如果蠕变指标为真，则 返回真值； 2. 否则dD=0。


FLAC3D采用面向对象的语言标准C++编写
本构模型都是以动态连接库文件(.DLL文件)的形式提供 VC++6.0(SP4)或更高版本的开发环境 优点
自定义的本构模型和软件自带的本构模型的执行效率
处在同一个水平
自定义本构模型(.DLL文件)适用于高版本的FLAC(2D)、
3DEC、UDEC等其他Itasca软件中
int iPlas = 0; double dTeTens = dTension; /* --- trial elastic stresses --- */ double dE11 = ps->stnE.d11; double dE22 = ps->stnE.d22; double dE33 = ps->stnE.d33;
初始化各式 中的常数
const char *UserMohrModel::Run(unsigned uDim,State *ps){ if ((uDim!=3)&&(uDim!=2)) return("Illegal dimension in Mohr constitutive model"); if(ps->dHystDampMult > 0.0) HDampInit(ps->dHystDampMult); /* --- plasticity indicator: */ /* store 'now' info. as 'past' and turn 'now' info off ---*/ if (ps->mState & mShearNow) ps->mState = (unsigned long)(ps->mState | mShearPast); ps->mState = (unsigned long)(ps->mState & ~mShearNow); if (ps->mState & mTensionNow) ps->mState = (unsigned long)(ps->mState | mTensionPast); ps->mState = (unsigned long)(ps->mState & ~mTensionNow);