【清华】实验6非线性方程求解

合集下载

实验六_非线性方程求解

实验六_非线性方程求解

1 ai ln P 0, i 1, 2,3, 4
求解上述非线性方程组即可得出均相共沸混合物的沸点和组成。 编写 M 文件 azeofun.m function f=azeofun(XT,n,P,a,b,c,Q) x(n)=1; for i=1:n-1 x(i)=XT(i); x(n)=x(n)-x(i); end T=XT(n); p=log(P); for i=1:n d(i)=x*Q(i,1:n)'; dd(i)=x(i)/d(i); end for i=1:n f(i)=x(i)*(b(i)/(T+c(i))+log(x*Q(i,1:n)')+dd*Q(1:n,i)-a(i)-1-p); end 编写程序生成 ai , bi , ci 和 Q n=4; P=760; a=[18.607,15.841,20.443,19.293]'; b=[2643.31,2755.64,4628.96,4117.07]'; c=[239.73,219.16,252.64,227.44]'; Q=[1.0,0.192,2.169,1.611;0.316,1.0,0.477,0.524;0.377,0.360,1.0,0.296;0.524,0.282,2.06 5,1.0]; 利用不同初值计算(仅展示有限结果) 编写程序 XT0=[0.5,0.5,0,90]; XT=fsolve(@azeofun,XT0,[],n,P,a,b,c,Q) 得到如下结果 XT =
在压强 p 760 mmHg 下,为了形成均相共沸混合物,温度和组分分别是多少?请 尽量找出所有可能的解。 共沸混合物,是指由两种或两种以上物质组成的液体混合物,当在某种压力下被蒸 馏或局部汽化时,在气体状态下和液体状态下保持相同的组分。 该混合物由 4 个组分组成,组分 i 所占的比例为 xi i 1, 2,3, 4 ,则

非线性方程求解实验

非线性方程求解实验

第二部分:提高实验(例题讲解) 第二部分:提高实验(例题讲解) 1、编写二分法法的Matlab程序。 2、利用所编写的程序求解非线性方程的近似解。 例如:求方程 f ( x ) = x 3 − x − 1 = 0在区间[1,1.5]内 的根,要求误差不超过0.005. 首先定义函数erfenfa.m 然后在Matlab窗口下输入: >> erfenfa 输入函数f(x)='x^3-x-1' 输入区间=[1,1.5] 请输入误差=0.005
第二部分:提高实验(学生操作) 第二部分:提高实验(学生操作) 1、编写二分法法的Matlab程序。 2、利用所编写的程序求解非线性方程的近似解。 利用二分法求方程 x 3 − 2 x 2 − 4 x − 7 = 0在区间 [3,4]内的根,要求误差不超过0.0005.
第三部分:创新实验(学生操作) 第三部分:创新实验(学生操作) 1、编写Newton迭代法的Matlab程序。 2、利用源程序求解非线性方程的近似解。 利用牛顿法求方程 x 2 − 2 x − 3 = 0 的近似根,取初 始值为4,要求误差不超过0.0001.
实验二: 实验二:非线性方程典型解法的数值实验 数学实验室
本次实验的目的: 本次实验的目的: 1、掌握基本的绘图命令Plot,ezplot。 掌握基本的绘图命令Plot,ezplot。 Plot,ezplot 2、编写Matlab程序,利用二分法求非线性方程的 编写Matlab程序, Matlab程序 近似解。 近似解。 3、在Matlab软件中,用迭代法求解非线性方程的 Matlab软件中, 软件中 近似解。 近似解。 4、编写Matlab程序,利用牛顿法求非线性方程的近 编写Matlab程序, Matlab程序 似解。 似解。

实验六__非线性方程求根

实验六__非线性方程求根

实验六 非线性方程求根6.1实验目的① 掌握二分法、牛顿迭代法等常用的非线性方程迭代算法;② 培养编程与上机调试能力.6.2算法步骤:6.2.1二分法计算()0f x =的二分法如下:① 输入求根取间[,]a b 和误差控制量ε,定义函数()f x .如果 ()()0f a f b <,转②;否则退出选用其它求根方法② 当a b ε->时,计算中点()/2x a b =+以及()f x 的值;分情况处理()f x ε<:停止计算,*x x =,转④()()0f a f x <:修正区间[,][,]a x a b →()()0f x f b <:修正区间[,][,]x b a b → ③ *2a b x += ④ 输出近似根*x6.2.2牛顿迭代法给定初始值0x ,ε为根的容许误差,η为()f x 的容许误差,N 为迭代次数的容许值.① 如果'0()0f x =或迭代次数大于N ,则算法失败,结束;否则执行② ② 计算010'0()()f x x x f x =- ③ 若10x x ε-<或1()f x η<,则输出1x ,程序结束;否则执行④ ④ 令01x x =,转向①6.3实验内容:求方程22()320x f x x x e =-+-=在0.5附近的根.6.4实验要求:(1)选择一种计算机语言设计出二分法和牛顿法的程序,并且选择不同的初值,观察所需的迭代次数和迭代结果;(2)分析二分法和牛顿法在非线性方程求根中的优缺点和收敛速度二分法简单易行,但只有线性收敛,且仅限于求实根;牛顿法也是一种简单的迭代法,具有二阶收敛速度(在单根邻近处)的特点,但对初值的选择比较苛刻,否则可能不收敛.实验程序:二分法:clcx=erfen(0,1,0.00001)function m=erfen(a,b,e)m=0;if f(a)*f(b)<0while abs(a-b)>ex=(a+b)/2;if abs(f(x))<em=x;break;endif f(a)*f(x)<0b=x;endif f(x)*f(b)<0a=x;endm=(a+b)/2;endendfunction r=f(x)r=x^3-3*x^2+2-exp(x);结果:x =0.4265牛顿法:x=newton(0.5,0.01,0.01,20)function m=newton(x0,e,g,N)m=0;for i=1:Nif df(x0)==0break;endx1=x0-f(x0)/df(x0);if abs(x1-x0)<e|abs(f(x1))<g m=x1;break;elsex0=x1;endendfunction r=f(x)r=x^3-3*x^2+2-exp(x);function r=df(x)r=3*x^2-6*x-exp(x);结果:x =0.4265。

清华大学数学实验6 非线性方程求解

清华大学数学实验6 非线性方程求解

a
x0 b x
a
f (x)
x0
f (b) 0
b x x0 : (a,b)中点
f (x0) 0 a1 a, b1 x0 f (x0) 0 a1 x0, b1 b
(a,b) (a1,b1) (an,bn) 区间每次缩小一半,n足
够大时,可确定根的范围
不足 收敛速度较慢
非线性方程迭代法的基本思想
y P1 h1
O
P2
C(x)
r1
r2
1 Q 2
x
h2
P1h1
P2h2
h12 x2 3
h22 (s x)2 3
s
x
实例2 均相共沸混合物的组分
均相共沸混合物(homogeneous azeotrope) 是由两种或两种 以上物质组成的液体混合物,当在某种压力下被蒸馏或局部汽 化时,在气体状态下和在液体状态下保持相同的组分(比例)
k 1
(qij表示组分i与组分j的交互作用参数 ,可以通过实验近似得到)
实例2 均相共沸混合物的组分
P i Pi
ln
Pi
ai
T
bi
ci
bi
T ci
ln
n
x j qij
j 1
n
j1
ln i
x j q ji
n
n 1 ln x jqij
j1
1
ai
ln
P
0
n
j1
将原方程 f (x) 0 改写成容易迭代的形式x (x) , 选
合适的初值x 0
,
进行迭代: x (x )
k 1
k
(k 0,1,2, )
例1
f (x) x2 x 14 0 x (x)

非线性方程的求解方法

非线性方程的求解方法

非线性方程的求解方法一、引言在数学领域中,非线性方程是指未知量与其对自身的各次幂、指数以及任意函数相乘或相加得到的方程。

求解非线性方程是数学中一个重要而又具有挑战性的问题。

本文将介绍几种常见的非线性方程求解方法。

二、牛顿迭代法牛顿迭代法是一种经典的非线性方程求解方法,它利用方程的切线逼近根的位置。

设f(x)为非线性方程,在初始点x0附近取切线方程y=f'(x0)(x-x0)+f(x0),令切线方程的值为0,则可得到切线方程的解为x1=x0-f(x0)/f'(x0)。

重复这个过程直到满足精确度要求或迭代次数达到指定次数。

三、二分法二分法是一种简单而又直观的非线性方程求解方法。

它利用了连续函数的中间值定理,即若f(a)和f(b)异号,则方程f(x)=0在[a, b]之间必有根。

根据中值定理,我们可以取中点c=(a+b)/2,然后比较f(a)和f(c)的符号,若同号,则根必然在右半区间,否则在左半区间。

重复这个过程直到满足精确度要求或迭代次数达到指定次数。

四、割线法割线法是一种基于切线逼近的非线性方程求解方法,它与牛顿迭代法相似。

由于牛顿迭代法需要求解导数,而割线法不需要。

设f(x)为非线性方程,在两个初始点x0和x1附近取一条直线,该直线通过点(x0,f(x0))和(x1, f(x1)),它的方程为y=f(x0)+(f(x1)-f(x0))/(x1-x0)*(x-x0),令直线方程的值为0,则可得到直线方程的解为x2 = x1 - (f(x1)*(x1-x0))/(f(x1)-f(x0))重复这个过程直到满足精确度要求或迭代次数达到指定次数。

五、试位法试位法是一种迭代逼近的非线性方程求解方法。

它利用了函数值的变化率来逼近根的位置。

设f(x)为非线性方程,选取两个初始点x0和x1,然后计算f(x0)和f(x1)的乘积,如果结果为正,则根位于另一侧,否则根位于另一侧。

然后再选取一个新的点作为下一个迭代点,直到满足精确度要求或迭代次数达到指定次数。

非线性方程组求解

非线性方程组求解

非线性方程组求解非线性方程组在科学、经济等领域中应用广泛,然而,由于非线性方程组的求解困难性,这使得许多问题存在困扰。

非线性方程组求解是一个复杂的过程,在此过程中需要对多种数学技术和算法有深入的了解。

本文就非线性方程组求解这个话题进行了探讨。

一、非线性方程组的定义非线性方程组是指一组包含至少一个非线性方程的方程组。

非线性方程组是一种数据的数学模型,它描述了在特定条件下各个因素之间的相互依赖关系。

非线性方程组的解通常用来预测一个系统的行为,并且是许多数学和科学领域的重要工具。

二、非线性方程组求解的困难性非线性方程组求解的困难性是因为它们存在着多个未知数和多个方程之间的相互依赖关系。

这使得非线性方程组的求解无法通过简单的代数运算来获得,而且通常需要更高级的数学知识和算法。

在许多情况下,非线性方程组可能无法解析地求解,这时需要采用数值方法来求解。

三、非线性方程组求解的方法1. 牛顿迭代法牛顿迭代法是最常用的求解非线性方程组的方法之一。

它将非线性方程组看作一组关于未知量的函数,并利用泰勒公式将其逼近为线性表达式。

由于直接求解非线性方程组比较难,牛顿迭代法通常将其转化为求解一系列线性方程组的问题。

2. 非线性迭代法非线性迭代法是一种通过递推计算的方式求解非线性方程组的方法。

具体地说,非线性迭代法会将非线性方程组转化为一组迭代公式,然后通过不断迭代来逼近方程组的解。

3. 二分法二分法是一种通过对非线性方程组的解进行区间逼近来求解的方法。

二分法的基本思路是通过每次将原来的区间对半分来寻找解所在的范围。

四、结语非线性方程组求解是一个重要的数学问题,应用广泛且具有挑战性。

本文主要介绍了三种很常用的求解方法,即牛顿迭代法、非线性迭代法和二分法。

在实际运用中,这些方法可以单独或者联合使用,以求得更准确的解。

非线性方程求解

非线性方程求解

非线性方程求解在数学中,非线性方程是一种函数关系,其表达式不能通过一次函数处理得到。

与线性方程不同,非线性方程的解决方案往往更具挑战性,因为它涉及到更复杂的计算过程。

尤其在实际应用中,非线性方程的求解是一个非常重要的问题。

本文将讨论几种常用的非线性方程求解方法。

二分法二分法,也称为折半法,是一种基本的求解非线性方程的方法之一。

它的核心思想是将区间一分为二并判断方程在哪一半具有根。

不断这样做直到最终解得精度足够高为止。

下面是利用二分法求解非线性方程的流程:1. 设定精度值和区间范围2. 取区间的中点并计算函数值3. 如果函数值为0或函数值在给定精度范围内,返回中点值作为精确解4. 如果函数值不为0,则判断函数值的正负性并缩小区间范围5. 重复步骤2-4直到满足给定精度为止当然,这种方法并不总是能够找到方程的解。

在方程存在多个解或者区间范围不合适的情况下,二分法可能会导致求解失败。

但它是一种很好的起点,同时也是更复杂的求解方法中的一个重要组成部分。

牛顿迭代法牛顿迭代法是一种更复杂的求解非线性方程的方法。

它利用泰勒级数和牛顿迭代公式,通过不断迭代来逼近根的位置。

下面是利用牛顿迭代法求解非线性方程的流程:1. 先取一个近似值并计算函数值2. 求出函数的导数3. 利用牛顿迭代公式,计算下一个近似根4. 检查下一个近似根的精度是否满足条件,如果满足,返回当前近似根5. 如果精度不满足,则将新的近似根带入公式,重复步骤2-5当然,牛顿迭代法的收敛性并不总是保证的。

如果迭代过程太过温和,它可能无法收敛到精确解。

如果迭代过程过于暴力,则会出现发散现象,使得求解变得不可能。

其他方法此外,还有一些其他的求解非线性方程的方法,例如黄金分割法、逆二次插值法、牛顿切线法等等。

其中每一种方法都有其优缺点,不同的情况下,不同的方法都可能比其他方法更加适合。

结论总体来说,求解非线性方程的方法非常复杂。

无论是哪种方法,都需要一定的数学基础和计算机知识。

【清华】实验6非线性方程求解

【清华】实验6非线性方程求解

实验六非线性方程求解姓名:李严凯学号:2009011887 班级:化93实验目的1掌握用MATLAB 软件求解非线性方程和方程组的基本用法,并对结果作出初步的分析2 练习用非线性方程组建立实际问题的模型并进行求解实验内容3.问题复述:(1)小张夫妇以按揭方式贷款买了一套价值20万的房子,首付了5万元,每月还款1000元,15年还清。

问贷款利率是多少?(2)某人欲贷款50万元购房,他咨询了两家银行,第一家银行开出的条件是每月还4500元,15年还清;第二家银行开出的条件是每年还45000元,20年还清。

从利率方面看,那家银行较优惠(简单假设年利率=月利率*12)?解答:模型:在实际问题中,购房人在支付首付后,向银行按揭申请贷款,并由银行支付购房人不足以支付的部分房价。

而后购房人将在一定时期内定期定额缴纳月供,直至还清银行的贷款本金及利息。

而银行则每月设定固定的月利率,按照复利的方式计算购房人应需缴纳的房款。

不妨设购房之初,购房人向银行申请了a0的贷款,银行的月利率为x,每月购房人向银行支付的月供为b,则在以后各月(第1、2、3…n个月),购房人所欠银行的贷款为:a1=a0*(1+x)-b;a2=a1*(1+x)-b;a3=a2*(1+x)-b;……an=a(n-1)*(1+x)-b当购房人所欠银行的贷款数an=0时,其贷款还清,还贷过程结束。

到此模型分析完毕,可以进行MATLAB求解。

求解:首先写出函数,表示出第n个月购房人所欠银行贷款:function y = interest( x,n ) %x为月利率n=180; %n为还清贷款的月数a(1)=15*(1+x)-0.1; %第1个月欠款for i=2:na(i)=a(i-1)*(1+x)-0.1; %第n个月欠款endy=a(n);end得到以上表达式后,令其中y=0即可解得月利率x,编程如下:[x,fv,ef,out]=fzero(@interest,0.05)考虑到银行贷款利率的实际情况,初值设为了0.05程序输出结果如下:x = 0.002081163889460fv = -5.287437154777308e-014ef =1out = intervaliterations: 12iterations: 14funcCount: 38algorithm: 'bisection, interpolation'message: 'Zero found in the interval [-0.014, 0.0952548]' 当然,也可以采用fsolve求解,列出:[x,fv,ef,out]=fsolve(@interest,0.05)得到x =0.002081163889460fv =-5.287437154777308e-014ef =1out = iterations: 14funcCount: 30algorithm: 'trust-region dogleg'firstorderopt: 9.769300386433235e-011message: [1x695 char]可以看出,两种解法的结果基本相同,不过fzero的解答时间更快,效果更优。

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

实验六非线性方程求解姓名:李严凯学号:2009011887 班级:化93实验目的1掌握用MATLAB 软件求解非线性方程和方程组的基本用法,并对结果作出初步的分析2 练习用非线性方程组建立实际问题的模型并进行求解实验内容3.问题复述:(1)小张夫妇以按揭方式贷款买了一套价值20万的房子,首付了5万元,每月还款1000元,15年还清。

问贷款利率是多少?(2)某人欲贷款50万元购房,他咨询了两家银行,第一家银行开出的条件是每月还4500元,15年还清;第二家银行开出的条件是每年还45000元,20年还清。

从利率方面看,那家银行较优惠(简单假设年利率=月利率*12)?解答:模型:在实际问题中,购房人在支付首付后,向银行按揭申请贷款,并由银行支付购房人不足以支付的部分房价。

而后购房人将在一定时期内定期定额缴纳月供,直至还清银行的贷款本金及利息。

而银行则每月设定固定的月利率,按照复利的方式计算购房人应需缴纳的房款。

不妨设购房之初,购房人向银行申请了a0的贷款,银行的月利率为x,每月购房人向银行支付的月供为b,则在以后各月(第1、2、3…n个月),购房人所欠银行的贷款为:a1=a0*(1+x)-b;a2=a1*(1+x)-b;a3=a2*(1+x)-b;……an=a(n-1)*(1+x)-b当购房人所欠银行的贷款数an=0时,其贷款还清,还贷过程结束。

到此模型分析完毕,可以进行MATLAB求解。

求解:首先写出函数,表示出第n个月购房人所欠银行贷款:function y = interest( x,n ) %x为月利率n=180; %n为还清贷款的月数a(1)=15*(1+x)-0.1; %第1个月欠款for i=2:na(i)=a(i-1)*(1+x)-0.1; %第n个月欠款endy=a(n);end得到以上表达式后,令其中y=0即可解得月利率x,编程如下:[x,fv,ef,out]=fzero(@interest,0.05)考虑到银行贷款利率的实际情况,初值设为了0.05程序输出结果如下:x = 0.002081163889460fv = -5.287437154777308e-014ef =1out = intervaliterations: 12iterations: 14funcCount: 38algorithm: 'bisection, interpolation'message: 'Zero found in the interval [-0.014, 0.0952548]' 当然,也可以采用fsolve求解,列出:[x,fv,ef,out]=fsolve(@interest,0.05)得到x =0.002081163889460fv =-5.287437154777308e-014ef =1out = iterations: 14funcCount: 30algorithm: 'trust-region dogleg'firstorderopt: 9.769300386433235e-011message: [1x695 char]可以看出,两种解法的结果基本相同,不过fzero的解答时间更快,效果更优。

由以上结果可以回答问题:银行的贷款月利率为0.21%。

(2)在本问中假定有年利率=月利率*12,因而比较两者的利率优惠程度,只需比较二者的月利率的大小。

模型分析的过程与第(1)问相同,只需将interest函数稍加修改即可。

对于第一家银行,将函数修改后为:function y = interest( x,n ) %x为月利率n=180; %n为还清贷款的月数a(1)=50*(1+x)-0.45; %第1个月欠款for i=2:na(i)=a(i-1)*(1+x)-0.45; %第n个月欠款endy=a(n);end其中还清贷款的月份仍为180个月,但每月还款变为了0.45万元。

运行fzero,有[x,fv,ef,out]=fzero(@interest,0.05)其中初值仍设为0.05运行结果为:x =0.005850792582847fv =1.616046185759501e-011ef = 1out = iterations: 13funcCount: 28algorithm: 'trust-region dogleg'firstorderopt: 1.722564940850133e-007message: [1x695 char]即第一家银行的月利率为0.585%对第二家银行,函数修改为:function y = interest( x,n ) %x为年利率n=20; %n为还清贷款的年数a(1)=50*(1+x)-4.5; %第1个年欠款for i=2:na(i)=a(i-1)*(1+x)-4.5; %第n个年欠款endy=a(n);end注意此时所求为年利率输入[x,fv]=fsolve(@interest,0.05)可得x =0.063948777092671fv =3.929887526510356e-010则按照题中假设,第二家银行的月利率应为x/12=0.0053290647577226,即0.532%故从月利率比较,第二家银行的利率小于第一家银行,则第二家银行更优惠。

易求得,按照第一家银行的规定,15年共需还款81万元,而第二家20年需还款90万元,即若按还款总额计算第二家要优于第一家银行。

另外,考查题中的假设条件“年利率=月利率*12”:若假定第二家银行实行月供,每月需还款3750元。

修改函数程序如下:function y = interest( x,n ) %x为月利率n=240; %n为还清贷款的月数a(1)=50*(1+x)-0.375; %第1个月欠款for i=2:na(i)=a(i-1)*(1+x)-0.375; %第n个月欠款endy=a(n);end利用fzero可求得其月利率为x =0.005479234984688,即为0.548%,与已求得的值0.532%差别不大;同样可求第一家银行的平均年利率为x =0.067402232214798,即6.74%;利用月利率*12求得年利率为7.02%,差别稍大。

故题中对于“年利率=月利率*12”的假设略显不足。

当然,如果只做定性比较,这样的精度可以容忍。

5.问题复述:由汽缸控制关闭的门,关闭状态的示意图如图1。

门宽a,门枢在H处,与H相距为b处有一门销,通过活塞与圆柱形的汽缸相连,活塞半径为r,汽缸长,汽缸内气体的压强为。

当用力F推门,使门打开一个角度α时(示意图如图2),活塞下降的距离为c,门销与H的水平距离b保持不变,于是汽缸内的气体被压缩,对活塞的压强增加。

已知在绝热条件下,气体的压强p和体积V满足p=c,其中γ是绝热系数,c是常数。

试利用开门力矩和作用在活塞上的力矩相平衡的关系(对门枢而言),求在一定的力F作用下,门打开的角度。

设a=0.8m,b=0.25m,r=0.04m,0.5m,=N/,γ=1.4,F=25N.解答:模型:本题中,基本的物理模型为杠杆平衡,同时涉及热力学的相关内容。

对于最终的平衡状态,设活塞向内移动了x m,应有:Fa cosα =(p-)sb ①其中,s=π同时有热力学的条件可知,=②其中,=s,V=s(-x),x/b=tanα解①②方程,代入相关关系式得Fa cosα =[-1]πb在不考虑大气压的情况下,原方程简化为Fa cosα =πb以上方程可由MATLAB编程求解求解:①先考虑简化情况,即不考虑大气压的影响:写出待求解函数:function y = balance1(al)y=10000*((0.5/(0.5-0.25*tan(al)))^1.4)*pi*0.04^2*0.25-25*0.8*cos(al); end以上是将各项数值代入的结果。

之后求解该方程。

调用fzero,有[x,fv]=fzero(@balance1,pi/6)其初值设为了30°结果为x =0.432945646735310fv =0则门打开的角度为α=0.4329rad,即α=0.4329*180/3.14°=24.8041°②再考虑较复杂的情况,即存在大气压的影响:同理,由以上过程可以先写出待求解函数:function y = balance1(al)y=10000*((0.5/(0.5-0.25*tan(al)))^1.4-1)*pi*0.04^2*0.25-25*0.8*cos(al );end以上是将各项数值代入的结果。

之后求解该方程。

调用fzero,有[x,fv]=fzero(@balance1,pi/6)其初值设为了30°结果为x=0.712003217036232fv=1.776356839400251e-015则门打开的角度为α=0.7120rad,即α=0.7120*180/3.14°=40.7958°7.问题复述:用迭代公式计算序列{},分析其收敛性,其中a分别取5,11,15;b(>0)任意,初值。

观察是否有混沌现象出现,并找出前几个分岔点,观察分岔点的极限趋势是否符合Feigenbaum常数揭示的规律。

解答:运用图形,分析在参数a取不同值时序列{}的收敛性计算的程序如下:a=[5,11,15];x=1; %赋初值n=50;for i=1:3A=a(i); %取a值for j=1:nx(j+1)=A*x(j)*exp(-x(j)); %迭代计算endxx( : ,i)=x';endk=(0:50)';[k,xx] %输出结果subplot(2,2,1),plot(k,xx(:,1)),title('a=5'), %画出3张图subplot(2,2,2),plot(k,xx(:,2)),title('a=11'),subplot(2,2,3),plot(k,xx(:,3)),title('a=15')得到的3张图如下:可以看出,在a=5时原数列收敛于(约)1.6,而在a=11,a=15时,原数列均不收敛。

且a=11时,原数列有两个收敛的子序列;a=15时,原数列有若干收敛的子序列。

以下考察混沌现象的出现:首先编写以下函数:function chaos( iter_fun,x0,r,n )kr=0;for rr=r(1):r(3):r(2)kr=kr+1;y(kr,1)=feval(iter_fun,x0,rr);for i=2:n(2)y(kr,i)=feval(iter_fun,y(kr,i-1),rr);endendplot([r(1):r(3):r(2)],y(:,n(1)+1:n(2)),'k.');end而迭代函数如下:function y = iter(x,a)y=a*x*exp(-x);020406011.52a=502040600246a=1102040600246a=15end输入命令如下:chaos(@iter,1,[5,20,0.01],[100,200])可得图像如下:5101520出现了混沌现象,读出前5个分岔点为,,,,,初步计算发现其比值分别为2.9,4.1,4.3,逐渐趋于Feigenbaum常数4.6692。

相关文档
最新文档