数学建模 第一章.docx
《数学建模 第一章.docx》由会员分享,可在线阅读,更多相关《数学建模 第一章.docx(23页珍藏版)》请在冰豆网上搜索。
![数学建模 第一章.docx](https://file1.bdocx.com/fileroot1/2022-11/15/3b519ebd-9742-41cb-9286-945e115bb1c2/3b519ebd-9742-41cb-9286-945e115bb1c21.gif)
数学建模第一章
第一章方程(组)模型
本章学习目的:
1.复习求解方程的基本原理和方法,掌握解方程的图形放大法和迭代算法;
2.能利用MATLAB软件编写迭代算法程序,了解迭代过程的图形表示;
3.熟练掌握用MATLAB软件的函数来求解方程和方程组;
4.通过范例展现求解实际问题的初步建模过程和MATLAB程序设计。
1.1引言
“方程是很多工程和科学工作的发动机”。
研究大型的土建结构、机械结构、输电网络、管道网络,研究经济规划、人口增长、种群繁殖等问题时,简单的分析可以直接归结为线性或非线性方程组,复杂一些要用到(偏)微分方程,求数值解时将转化为n非常大的方程组。
若干世纪以来,工程师和科学家花了大量的时间用于求解方程(组),数学家研究各种各样的方程求解方法。
本章我们就是要学习求解线性方程组、非线性方程(组)的方法,以及利用数学软件利用计算机对方程和方程组进行求解。
1.2方程的求解方法
考虑求方程f(x)=0的解,我们通常采用这样的几种方法:
因式分解法、图形放大法、数值迭代逼近法。
1.因式分解法:
这是我们最熟悉、常用的一种方法,这个方法的关键在分解因式,包括对多项式函数、三角函数和指数函数等的分解。
但对于无法进行分解的函数则无能为力。
2.图形放大法:
由于计算机的广泛应用,可以非常方便地作出函数f(x)的图形(曲线),找出曲线与x轴的交点的横坐标值,就可求出f(x)=0的近似根。
这些值尽管不精确,但是直观,方程有多少个根、在什么范围,一目了然。
并且可以借助于计算机使用图形局部放大功能,将根定位得更加准确。
3.数值迭代逼近法:
利用图形的方法或连续函数的零点存在性定理,可以推知f(x)在某一区间内有根,我们就可以用数值方法来求方程的根,这就是迭代逼近法。
迭代逼近法分为区间的迭代和点的迭代。
区间迭代又分为对分法和黄金分割法;点的迭代又分为简单迭代法、单点割线法、两点割线法、牛顿法等。
迭代失败后又可以采用加速迭代收敛方法。
1.2.1图形放大法
用图形放大法求解方程f(x)=0的步骤:
(1)建立坐标系,作曲线f(x);
(2)观察f(x)与x轴的交点;
(3)将其中的一个交点进行局部放大;
(4)该交点的横坐标值就是方程的一个根;
(5)对所有的交点进行相同的处理,就得到方程的所有解。
例1.1求方程
所有的根及大致分布范围,欲寻求其中的一个实根,并且达到一定的精度。
(1)画出
的图形;
x=-6:
0.01:
6;
y=x.^5+2*x.^2+4;
plot(x,y)
gridon;%画坐标格
我们可以看出方程在-2~+2范围有一个实根。
(2)逐次缩小范围得到较精确的根。
x=-2:
0.01:
2;
y=x.^5+2.*x.^2+4;
plot(x,y)
gridon
x=-2:
0.01:
-1.5;
y=x.^5+2.*x.^2+4;
plot(x,y)
gridon
x=-1.6:
0.01:
-1.5;
y=x.^5+2.*x.^2+4;
plot(x,y)
gridon
因此我们可以看出这个实根的值在-1.55~-1.5之间。
1.2.2简单迭代法
1.迭代算法步骤:
对方程f(x)=0求解
(1)对方程经过简单变形得到
(不是唯一的),x被称之为不动点;
(2)设置
为迭代初值,迭代过程为
,n=0,1,2……
(3)当两次迭代结果之差小于某个设定的误差值时,我们认为迭代结果是收敛的,可得到结果的近似值
。
例1.2求方程
的非负实根。
解:
由于函数
连续,并且在x=0和x=1处函数值符号相反,可以判断函数在区间(0,1)必有零点,即方程
在(0,1)内必然存在根。
(1)先将函数变形为
;
(2)设置迭代初值为0,编程进行迭代。
n=1;
x=0;
y=exp(x)/3;
ys=vpa(y,10);%给出y的数值型结果,有效位数为10
z=abs(y-x);
whilez>10^(-5)
x=y;
y=exp(x)/3;
ys=vpa(y,10);
z=abs(y-x);
n=n+1;
end
n,y,ys
n=
21
y=
0.6190
ys=
.6190471917
从该结果可以看出,迭代21次后两次迭代的结果误差值满足小于
的条件,结果收敛,迭代结果为0.6190,若保留小数点后10位有效数字则结果为0.6190471917。
例1.3用迭代方法求解方程
解:
(1)对方程变形为
,有不同的形式,比如
;(a)
;(b)
;……(c)
(2)设定初始值为1,编程迭代求解
x=1;y=1;z=1;
fork=1:
25
x=x^3-x^2-1;
y=(y^2+y+1)^(1/3);
z=1+1/z+1/z^2;
end
x,y,z
x=
-Inf
y=
1.8393
z=
1.8393
在程序中,函数x,y,z分别对应方程(a)(b)(c),从结果可以看出方程(a)不收敛,结果趋于负无穷大,方程(b)(c)收敛,结果为1.8393。
而且,还可证明(b)比(c)收敛速度快。
注:
这段程序和例1.2有所不同,这里是设定了固定的迭代次数。
2.迭代失败的改进:
加速迭代收敛方法
例1.3中方程(a)的迭代是失败的(即迭代不收敛),如何解决?
当我们遇到迭代失败时,可以采用以下的简单方法解决。
我们考虑不直接用
迭代,而用
和x的加权平均
进行迭代,
为参数。
即:
在满足
的条件下,取
,由此可以得到
。
在实际迭代过程中用
代替a,因此在迭代中
,
其加速迭代过程如下:
。
采用以上方法对例2.3中方程(a)进行改进,产生新的迭代函数
加速迭代过程为
。
再编写程序计算:
x=0;
fork=1:
20
x=(-2*x^3+x^2-1)/(-3*x^2+2*x+1);
end
x
x=
1.8393
x=100;
fork=1:
20
x=(-2*x^3+x^2-1)/(-3*x^2+2*x+1);
end
x
x=
1.8393
由此我们看到选用两个不同的初值0和100都能得到结果是1.8393,改进是非常有效的,同学们还可尝试不同的初值,观察其迭代的收敛性。
实验表明,它比(b)、(c)的收敛速度都快。
但要注意,x=1不能作为初值,这会出现被0除的错误。
通过这个例子我们看出,求解方程的迭代函数构造形式多样,不同迭代函数的收敛性和收敛速度都可能不同,需要在遇到具体问题时灵活应用。
1.3方程组的求解方法
1.3.1线性方程组的求解
我们在线性代数中已经学习了线性方程组的求解方法。
对于线性方程组
可以写成矩阵的形式
由线性代数的知识可知,线性方程组的解可能出现三种情形:
无解、有唯一解、有无穷多组解。
这主要取决于系数矩阵A的秩与增广矩阵(A︱b)的秩是否相等、秩与变量个数是否相等,具体地:
若R(A)≠R(A︱b),则
无解;
若R(A)=R(A︱b)=n(n为变量个数),则
有唯一一组解;
若R(A)=R(A︱b)有无穷多组解。
求矩阵A的秩可以很方便的用MATLAB的rank(A)函数求得。
求解线性方程组的方法大致可以划分为两类:
直接消去法、迭代数值解法。
直接消去法在线性代数中已经学过,这里不再赘述。
线性方程组可以看成是非线性方程组的特例,其迭代数值解法相同,在非线性方程组的迭代解法中介绍方法。
1.3.2非线性方程组的迭代解法
非线性方程组的一般形式为
可以改写为等价的方程组
用这个方程组进行迭代求得精确解。
例1.4求解方程组
解:
对方程组进行变形,构造如下的迭代函数:
(1)
或
(2)
思考:
迭代序列如何表示?
求解方程组迭代产生的序列是数组:
对本例选择初始点(0,0),(2,3),(8,9),…分别计算,迭代次数逐渐增加,观察结果。
迭代程序如下:
(初始值是(2,3),迭代次数为20次)
x=[2,3];y=[2,3];
fork=1:
20
a=0.1*x
(1)^2+0.1*x
(2)^2+0.8;
x
(2)=0.1*x
(1)*x
(2)^2+0.1*x
(1)+0.8;
x
(1)=a;
b=(10*y
(2)-8)/(y
(2)^2+1);
y
(2)=sqrt(10*y
(1)-8-y
(1)^2);
y
(1)=b;
end
x,y
x=
1.00001.0000
y=
2.19343.0205
从这个结果看出,x
(1)=1,x
(2)=1,和y
(1)=2.1934,y
(2)=3.0205是方程组的两组解。
问:
程序中a、b变量的作用是什么?
取消可以吗?
1.4MATLAB软件直接求解法
1.4.1任意函数方程与线性方程组可用命令solve()求解
solve()语句的用法:
1.单变量方程:
f(x)=0
1)符号解:
例1.5求解方程:
解:
输入:
x=solve('a*x^2+b*x+c')
输出为:
x=
[1/2/a*(-b+(b^2-4*a*c)^(1/2))]
[1/2/a*(-b-(b^2-4*a*c)^(1/2))]
或输入:
solve('a*x^2+b*x+c')
输出为:
ans=
1/2/a*(-b+(b^2-4*a*c)^(1/2))
1/2/a*(-b-(b^2-4*a*c)^(1/2))
2)数字解:
如果不能求得精确的符号解,可以计算可变精度的数值解。
例1.6解方程:
解:
s=solve('x^3-2*x^2=x-1')
s=
1/6*(28+84*i*3^(1/2))^(1/3)+14/3/(28+84*i*3^(1/2))^(1/3)+2/3
-1/12*(28+84*i*3^(1/2))^(1/3)-7/3/(28+84*i*3^(1/2))^(1/3)+2/3+1/2*i*3^(1/2)*(1/6*(28+84*i*3^(1/2))^(1/3)-14/3/(28+84*i*3^(1/2))^(1/3))
-1/12*(28+84*i*3^(1/2))^(1/3)-7/3/(28+84*i*3^(1/2))^(1/3)+2/3-1/2*i*3^(1/2)*(1/6*(28+84*i*3^(1/2))^(1/3)-14/3/(28+84*i*3^(1/2))^(1/3))
double(s)
ans=
2.2470
-0.8019+0.0000i
0.5550-0.0000i
vpa(s,10)
ans=
2.246979604+.1e-9*i
-.8019377358-.1866025404e-9*i
.5549581322-.1339745960e-10*i
该方程无实根。
3)无穷解
例1.7求解方程:
解:
输入:
solve('tan(x)-sin(x)=0')
输出为:
ans=
0
该方程有无穷多个解,不能给出全部解,这里只得到其中的一个。
2.多变量方程组:
例1.8解方程组
解:
输入:
[x,y]=solve('x^2*y^2','x-y/2-b')
x=
0
0
b
b
y=
-2*