MATLAB论文12010245327马文建.doc

《MATLAB》课程论文 利用MATLAB语言解决力学中非线性问题 姓 名：马文健 学 号：12010245327 专 业：电子信息工程 班 级：2010级电子班 指导老师：汤全武 学 院：物电学院 完成日期：2011-12-14 利用Matlab语言解决力学中非线性问题 （马文健 12010245327 2010级电子班） [摘要]Matlab 因其强大的功能而在世界范围内得以流行，特别在工程领域的应用甚为广泛。本文从弹性力学中的一个问题出发，先简要叙述并分析了目前厚壁筒弹塑性应力分析中所遇到的问题及其Matlab 在其中应用的思考；后建立了厚壁筒弹塑性区的模型，进行了各应力分析，用Matlab 中介绍的相关方法求解了弹塑性分界面半径并绘制了各应力关于半径的变化曲线，最后根据模型求解的结果进行了简要的分析。 [关键词]Matlab 弹塑性应力 厚壁筒 一、提出问题 厚壁筒在承受内压载荷的作用下，随着压力的增加，筒壁应力不断增加。厚壁筒在承受逐渐增加压力的过程中，会经历弹性阶段、筒体部分屈服阶段、整体屈服阶段、材料硬化、筒体过度变形、直至爆破失效阶段。而在分析厚壁筒的弹塑性应力分布时，遇到两个问题： 【1】弹塑性区分界面的半径的精确确定较为困难，因此半径是一个非线性方程； 【2】弹塑性应力的分布曲线绘制难以精确，不能真实反映不同半径处的应力状况。 二、分析问题 针对厚壁筒弹塑性应力分析时遇到的两个问题，可利用Matlab 软件中的相关知识解决。对于弹塑性区分界面半径的确定，可利用方程求根的相关方法，从二分法、开方法、Newton 法、Newton 下山法以及弦截法中选取精度和收敛速度均较佳的方程求根数值方法进行求解；同时还可以利用Matlab 中的Solve 函数求其精确解。 对于弹塑性应力分布曲线的绘制，可以根据求出的各应力分布的具体方程，用Plot 函数进行绘制，同时要注意曲线后期的处理工作，以便更好地从图中得到各应力随半径分布的情况。 1、模型建立及求解 图1 受内压厚壁圆筒 图2 弹性区 图3 塑性区 为简化分析，假设厚壁圆筒为理想弹塑性体，不考虑材料在塑性变形过程中塑性强化，筒体仅受内压 pi 作用，筒体的内半径为Ri ，外半径为Ro 。 初始假设厚壁圆筒的内半径为Ri ?152.5mm，外半径为Ro ??254mm，内壁所受的压力为pi=340MPa，取其屈服极限为=750MPa 先分析筒体仅受内压Pi作用的情形。当内压pi 大于弹性极限压力pe 时，圆筒内壁的屈服区向外扩展，筒体沿壁可以分成塑性区和弹性区两个区域，其中内侧为塑性区，外侧为弹性区。弹性区的内外半径分别为Rc 和Ro ，承受内压力为Pc ，见图2；假想两区域的交界圆面的半径为Rc ，则塑性区的内外半径分别为Ri 和Rc ，承受的内外压力分别为Pi和Pc ，如图3所示。 1.1 塑性区（RirRc ） 材料处于塑性状态时，筒壁微元体的平衡微分方程成立，即。 假设材料塑性变形时符合Tresca屈服条件，代入上式，则有：，积分上是既有。易知有边界条件：。将第一个条件代入式[4]。可求出A，并代入Tresca屈服条件，可得到塑性区各应力分量的表达式： 将第二个边界条件代入式[7]中可得塑性区交界面压力为。 1.2 弹性区（coR . r . R） 1.2.1应力分析 弹性区内壁面即为弹塑性区交界面，即弹性区内壁面呈塑性状态。由式[7]、[8]、[9] 可知在r=Rc有弹性区内壁面处各应力表达式为: 弹性区内壁面开始屈服时，假设应力符合Tresca屈服条件，则有。 由连续性可得式[10]和式[4]应为同一值，即有： 1.2.2 计算流程 其中??为下山因子，??为精度指标， 0 x 为迭代初值， N 为预设迭代次数。 ???0 ????1?的选择是个逐步探索的过程，从???1开始反复将因子??的值减半进行计算，一旦单调性条件式[22]成立，则下山成功。 ??在此可取为???10 ^ (?3)即可达到要求。 1.2.3 程序编写 （1）编写文件名为Mendnewton.m 的Matlab 文件 function [x,k]=Mendnewton(f,x0,emg) %输入形参表 [f1,d1]=feval(f,x0); k=1; %迭代次数初值 x(1)=x0; %迭代初值赋值于 x(2)=x(1)-f1/d1;