【2019年整理】岩土工程数值计算.ppt

【2019年整理】岩土工程数值计算.ppt

《【2019年整理】岩土工程数值计算.ppt》由会员分享，可在线阅读，更多相关《【2019年整理】岩土工程数值计算.ppt（227页珍藏版）》请在冰豆网上搜索。

岩土工程数值计算二零一一年三月参考文献参考文献计算土力学上海科学技术出版社朱百里土工计算机分析中国建筑工业出版社龚晓南土工数值计算中国铁道出版社钱家欢土工原理与计算水利电力出版社钱家欢岩土工程有限元分析理论与应用科学出版社谢康和FiniteelementanalysisingeotechnicalengineeringDavidM.Potts有限元法的基本知识有限元法的基本知识0.11自由度的编号自由度的编号N=n+1READ（1,*）（IB3（I）,I=1,NP3）READ（1,*）（LH（I）,I=1,NP5）DO60I=1,NP3IS=IB3（I）CALLHDIR（IS）DO60J=2,3IP3（JP

（1）,J-1）=JP（J）60CONTINUEDO70I=1,NP5IS=LH（I）/10II=LH（I）-IS*10IP3（IS,3）=II70CONTINUEN3=0DO80I=1,NPDO100J=1,3IF（IP3（I,J）.EQ.0）GOTO100IP3（I,J）=IP3（I,J）+N3N3=N3+1100CONTINUE80CONTINUE21301节点号X方向约束信息Y方向约束信息1210节点号0代表孔压22半带宽半带宽dmdmdm=b+1b为相邻自由度编号的最大差值（对应节点在同一单元）nnmmkkkkkkkkkkkkkkkkkK99887766555453524544353325221100dmIJK=3NN=12DO10I=1,4IVP=IVEN（I）DO10J=1,IJKIJ=IJK*（I-1）+JICN（IJ）=IP（IVP,J）10CONTINUEDO20J=1,NNICNP=ICN（J）DO20K=1,NNICNQ=ICN（K）IF（ICNP.EQ.0.OR.ICNQ.EQ.0）GOTO20IF（IA3（ICNP）.GT.（ICNP-ICNQ）GOTO20IA3（ICNP）=ICNP-ICNQ20CONTINUEHALFBANDWIDTH-1RELATEDTOIthFREEDEGREE33总刚的组装总刚的组装节点平衡法12345每一单元节点力与节点位移之间的关系为22222RFFFFe321333231232221131211321kkkkkkkkkFFF5214212213212211211kkkkkkF3212211213212211212kkkkkkF2213215213212211213kkkkkkF根据外荷载与节点力平衡条件直接刚度法先把每个单元的单刚阶数扩大成总刚阶数，把单刚中按局部编码的子块搬到总刚中相应的总体编码的位置中去，余下部分用零子块填充。

12345632100000000000000000000000000654321333231232221131211kkkkkkkkkK12345612344总刚的存储总刚的存储一维变带宽存储由于整体刚度矩阵具有对称性、稀疏性和带状性，采用变带宽下三角一维存储。

887877686766565546454435343324232214131211kkkkkkkkkkkkkkkkkkkkk887877686766565546454435343324232214131211kkkkkkkkkkkkkkkkkkkkkDO20I=1,NPDO20J=1,JANX=JA*（I-1）+JM=INE（I）20MA（NX）=IWU（M,J）M=NP*JADO60I=1,MK=MA（I）IF（K.LE.0）GOTO60DO40J=1,ML=MA（J）IF（L.EQ.0.OR.L.GT.K）GOTO40II=IDK（K）+L-KTK（II）=TK（II）+EK（I,J）40CONTINUE60CONTINUE将某一单元的Ke集整到劲度矩阵K中SKYLINE（轮廓线法）按列存储刚度矩阵上三角区必要部分，按行存储下三角区中必要部分。

对刚度矩阵中出现少数非常长的列的情况下，存储要求不会剧烈增加，很容易利用向量点积例行程序。

8878776867665856554846454438353433282423221814131211kkkkkkkkkkkkkkkkkkkkkkkkkkK11k12k22k13k23k33k14k24k34DO355K=1,4NRCC=3*（K-1）NR=NQ（NOD（K,IX）-1DO350M=1,2+IFLOWNRCC=NRCC+1NR=NR+1IF（ICODE.LT.2）THENLENNC=（LOCC（NR）-LOCC（NR-1）-1）/2DO345L=1,4NCCC=3*（L-1）NCN=NQ（NOD（L,IX）-1DO344N=1,2+IFLOWNCCC=NCCC+1NCN=NCN+1IF（NR.LT.NCN）NN=LOCC（NCN）-NCN+NRIF（NR.EQ.NCN）NN=LOCC（NR）IF（NR.GT.NCN）NN=LOCC（NR）-LENNC-NR+NCNS（NN-ISHIFT）=S（NN-ISHIFT）+C1（NRCC,NCCC）344CONTINUE345CONTINUEENDIFSL（NR）=SL（NR）+ZY（NRCC）350CONTINUE355CONTINUE（NCN,NR）（NR,NCN）55边界条件的引边界条件的引入入划0置1法：

处理ui=0约束乘大数法：

处理ui=R约束将总刚相应的主对角元素改为1，将对应的行、列其它元素改为0。

将荷载向量中相应的元素改为0。

或预先将每个节点的方程编号，已知位移的节点不编号。

将总刚相应的主对角元素置一大数，将荷载向量中相应的元素改为该大数乘R。

初等变换法：

处理两个变量有确定关系如轴对称问题，可以取一夹角为扇形区计算。

在斜对称AB边上的点只能沿AB上移动，所以点C上的点的位移有使用初等变换法消除一个相关方程，若消除对应的方程，可以将所在的行乘tg加到u所在的行，并将所在的列也乘tg加到u所在列上去。

tguccJP

（1）:

NODENUMBERJP

（2）:

X-DIRECTIONRESTRAININFORMATIONJP（3）:

Y-DIRECTIONRESTRAININFORMATIONJP

（1）=I3/100K3=I3-JP

（1）*100JP

（2）=K3/10K3=K3-JP

（2）*10JP（3）=K366有限元法解题步骤有限元法解题步骤1建立计算网格2设定计算相关参数（计算精度、控制等）3边界条件和内约束，力、位移、孔压等4选择单元类型、本构模型，并输入本构参数5计算结果提取、分析等岩土工程问题分析方法岩土工程问题分析方法1.1.11.1岩土工程问题控制方程的建立岩土工程问题控制方程的建立0uzyxzxyxx0vzyxzyyxy0gwzyxzyzxz11土体平衡方程土体平衡方程22土体本构方程土体本构方程Ddzzppdxxppgzxwwzpp33土体几何方程土体几何方程xuxyvyzwzxvyuxyywzvyzzuxwzx44土体有效应力原理土体有效应力原理pxxpyypzz55孔隙流体平衡方程孔隙流体平衡方程0ukvgxpwxxw0vkvgypwyyw0gwkvgzpwwzzw66渗流连续方程渗流连续方程zwyvxutzvyvxvzyxdzzzzdzxxxgzxw77总控制方程总控制方程uxpzxwddyxvddzudyudxud2551324412225522442211）（）（vypzywddyxuddzvdyvdxvd2551324412225522112244）（）（wgzpzyvddyzuddzwdywdxwd2551325513223322552255）（）（0wkzvkyukxzwyvxutggzpkzypkyxpkxzyxwwwzyx1.21.2岩土工程基本分析方法岩土工程基本分析方法11总应力分析法及其控制方程总应力分析法及其控制方程总应力分析法与一般固体力学相同。

从应用上讲，一般用于不考虑渗流固结的情况，如饱和粗粒土地基、透水性土料组成的土坝路堤的应力和变形分析以及饱和软粘土地基短期变形和稳定性分析。

0）（）（2551324412225522442211zxwddyxvddzudyudxud0）（）（2551324412225522112244zywddyxuddzvdyvdxvd0）（）（2551325513223322552255zyvddyzuddzwdywdxwd22有效应力分析法及其控制方程有效应力分析法及其控制方程在有效应力分析法中，土体的有效应力和孔压被严格区分，并将土骨架变形与孔隙水的渗透同步考虑。

因此，有效应力分析法较能更真实地反映土体的自身特性，能更合理地计算土体对载荷的响应，应用范围更广。

有效应力分析法尚需要有效应力原理和连续性方程。

总应力法中只有位移变量且仅与空间有关，而有效应力法中还有孔压变量，而且与空间和时间均有关。

0）（）（2551324412225522442211xpzxwddyxvddzudyudxud0）（）（2551324412225522112244ypzywddyxuddzvdyvdxvd0）（）（2551325513223322552255gzpzyvddyzuddzwdywdxwdzwyvxutggzpkypkxpkwwzyx22222233总应力分析法总应力分析法和和有效应力分析法有效应力分析法关系关系总应力法是有效应力法中当孔压p=0时的特殊形式。

在有效应力分析中如果采用与总应力分析相同的土工参数，并令孔压p=0，所得结果即为总应力分析结果。

总应力分析一般采用土体得不排水指标，由此进行的是加荷瞬时或短期应力和变形分析。

但也可采用土体的排水指标，此时进行的是最终或长期应力和变形分析。

当进行线弹性分析并采用排水指标，总应力分析得到的结果为有效应力分析的最终（孔压消散完毕、主固结完成）结果。

44不排水孔压计算不排水孔压计算f000ffffpppDffDfDDD333300IIKDef有效应力原理其中fsffpKnpKn1vv）（ezyxefKKpsfeKnKnK11孔压增量不仅使孔隙流体压缩，也会引起土颗粒的体积压缩。

相应的有效应力增量也会引起土颗粒的体积变化，然而，由于有效应力必须通过颗粒接触，但接触面积很小，导致体积变化也很小，如果忽略这种体积变化。

则总的体积变化为对于饱和土体，Ks、Kf都比土骨架模量大很多，在不考虑具体值时（此时准确值已不重要），可假设Kf=Ks,得feKKnKKfe由于Ks比土骨架模量大很多，如孔隙流体压缩性较大，以致Kskf，则对于排水分析，取Ke=0，加载过程中孔压不变。

）21（AAu3）1（）21（）1（A对于各向同性线弹性土体进行不排水分析时，必须设置Ke。

据经验，对于饱和土体，只要ke足够大，土体对ke的实际大小并不敏感。

但Ke太大时，可能导致数值不稳定，即不排水泊松比接近0.5。

专家建议设置ke=Kskel，在100和1000之间，Kskel是土骨架的体积模量。

101001000u=0.10.45200.49460.4994=0.30.47930.49770.4998）21（AAu3）1（）21（）1（A1.31.3岩土工程问题的边界条件岩土工程问题的边界条件11固结分析中的边界条件固结分析中的边界条件对节点位移和孔压已知的情况，可以有两种处理方法

（1）仍给以自由度编号，在解方程时把