上机作业文档格式.docx
《上机作业文档格式.docx》由会员分享，可在线阅读，更多相关《上机作业文档格式.docx（39页珍藏版）》请在冰点文库上搜索。
disp（'
不动点迭代法'
）;
n0=100;
p0=-1.5;
fori=1:
n0
p=（5^p0-4）/3;
ifabs（p-p0）<
=10^（-6）
ifp<
|p-p0|='
）
disp（abs（p-p0））
不动点迭代法求得方程的负根为:
'
disp（p）;
用不动点迭代法迭代步数为：
disp（i）;
break;
else
p0=p;
end
p1=1.5;
p=（log（3*p1+4））/log（5）;
ifabs（p-p1）<
ifp>
|p-p1|='
disp（abs（p-p1））
用不动点迭代法求得方程的正根为：
p1=p;
Newton迭代程序：
牛顿法'
n0=80;
p0=1.5;
p=p0-（3*p0-5^p0+4）/（3-log（5）*（5^p0））;
用牛顿法求得方程的正根为:
p1=-1.5;
p=p1-（3*p1-5^p1+4）/（3-log（5）*（5^p1））;
用牛顿法求得方程的负根为:
ifi==n0
disp（n0）
用牛顿迭代后无法求出方程的解'
运行结果：
结果分析：
选用不动点迭代公式x=（5^x-4）/3求负根的运行结果是收敛的，求正根发散，这是由于指数函数的存在使迭代发散。
而使用的迭代公式x=log（3x+4）/log5在求负根时发散，求正根时收敛，这是由于公式中存在对数函数。
由以上运算结果可知，Newton迭代法公式唯一，且求正负根均收敛，迭代步数明显少于不动点迭代，运算效率高。
2、用Newton法与重根计算法求解方程x-sinx＝0的根.再用Steffensen’smethod加速Newton法收敛，比较结果。
Newton迭代法原理：
重根计算法原理：
在Newton迭代法的基础上改进迭代公式，使得x（n+1）=x（n）－mf（x（n））/f'
（x（n）），m为根的重数。
牛顿法：
x=x-（x-sinx）/（1-cosx）,给定初值：
x=0
程序：
p0=1;
p=p0-（p0-sin（p0））/（1-cos（p0））;
用牛顿法求得方程的根为:
用Newton法迭代步数为：
重根计算法程序：
重根计算法'
m=2;
p=p0-m*（p0-sin（p0））/（1-cos（p0））;
用重根计算法求得方程的根为:
用重根计算法迭代步数为：
Steffensen加速Newton法程序：
Steffensen加速牛顿法'
i=0;
p
（1）=p0;
whilei<
=n0
fork=1:
2
p（k+1）=p（k）-（p（k）-sin（p（k）））/（1-cos（p（k）））;
p1=p
（1）-（p
（2）-p
（1））^2/（p（3）-2*p
（2）+p
（1））;
p=p1-（p1-sin（p1））/（1-cos（p1））;
用Steffensen加速牛顿法求得方程的根为:
用Steffensen加速Newton法迭代步数为：
i=i+1;
disp（i+1）;
p
（1）=p1;
由运行结果可知，Newton法收敛较慢，重根计算法效率大大提高，而Steffensen加速牛顿法效率最高，仅两次迭代就达到了要求的计算精度。
上机作业2：
分别用Jacobi、Seidel、Sor（w=0.8,1.1,1.2,1.3,1.4,1.5）方法求解方程组Ax=b。
这里
X0=[0……0]’;
e=10-5
要求：
抄题，公式，程序、计算结果（终止迭代步数k、近似解x（k）），结果分析（三种迭代列表）。
实验原理：
1、雅可比迭代法（J-迭代法）：
线性方程组
，可以转变为：
迭代公式
其中
，称
为求解
的雅可比迭代法的迭代矩阵。
以下给出雅可比迭代的分量计算公式，令
，由雅可比迭代公式有
，既有
，于是，解
的雅可比迭代法的计算公式为
2、高斯-赛德尔迭代法（GS-迭代法）：
GS-迭代法可以看作是雅可比迭代法的一种改进，给出了迭代公式：
其余部分与雅克比迭代类似。
3、逐次超松弛迭代法（SOR-迭代法）：
选取矩阵A的下三角矩阵分量并赋予参数w，将之作为分裂矩阵M，
，其中，w>
0，为可选择的松弛因子，又
（1）公式构造一个迭代法，其迭代矩阵为
从而得到解
的逐次超松弛迭代法。
其中：
由此，解
的SOR-迭代法的计算公式为
（2）
观察
（2）式，可得结论：
当w=1时，SOR-迭代法为J-迭代法。
当w>
1时，称为超松弛迭代法，当w<
1时，称为低松弛迭代法。
Jacobi迭代法程序：
用Jacobi迭代法求解线性方程组'
i=1;
i<
=100;
a=[-13111111111;
1-1311111111;
11-131111111;
111-13111111;
1111-1311111;
11111-131111;
111111-13111;
1111111-1311;
11111111-131;
111111111-13];
d=diag（diag（a））;
l=d-tril（a）;
u=d-triu（a）;
d0=inv（d）;
b=[-1;
-1;
-1];
x0=zeros（10,1）;
B=-d0*（l+u）;
f=d0*b;
x=B*x0+f;
whilenorm（x-x0,inf）>
=1e-5
x0=x;
x=B*x0+f;
i=i+1;
方程组的解为:
x
迭代次数:
i
Gauss-Seidel迭代法程序：
用Gauss-Seidel迭代法求解线性方程组'
B=-inv（d+l）*u;
f=inv（d+l）*b;
Sor迭代法程序：
分别令w=0.8,1.1,1.2,1.3,1.4,1.5
用Sor迭代法求解线性方程组'
w=0.8;
B=inv（d+w*l）*[（1-w）*d-w*u];
f=w*inv（d+w*l）*b;
松弛因子为:
w
运行结果分析：
由运行结果截屏的对比可发现，Jacobi迭代法收敛速度较慢，seidel
迭代法与松弛因子为1.1的Sor法的迭代效率一样，而w=0.8时，Sor的效率最高，随着松弛因子的增大，计算效率下降。
上机作业3
用Newton法与最速下降法求方程组x2+2y2-2=0,x2=y在（0.8，0.7）附近的根。
1）抄题；
2）迭代公式（初值）或简单原
理；
3）程序，结果；
（打印）4）结果分析。
牛顿迭代法原理
对于非线性方程
在x（k）处按照多元函数的泰勒展开，并取线性项得到
初值取x=0.8,y=0.7;
实验步骤
步骤一：
编写牛顿迭代法的基本程序。
打开Editor编辑器，输入以下语句：
function[x,n,data]=new_ton（x0）
ifnargin==1
tol=1e-6;
x1=x0-f1（x0）/df1（x0）;
n=1;
while（norm（x1-x0）>
tol）
x0=x1;
n=n+1;
data（:
n）=x1;
x=x1;
以文件名new_ton.m保存。
步骤二：
编写方程函数与方程的Jacobi矩阵函数。
（1）打开Editor编辑器输入以下语句：
%牛顿迭代法的方程函数
functionf=f1（x0）
x=x0
（1）;
y=x0
（2）;
f1=x^2-2*x-y+0.5;
f2=x^2+4*y^2-4;
%最后方程函数以行向量输出
f=[f1f2];
以文件名f1.m保存。
（2）新打开Editor编辑器输入以下语句：
%牛顿迭代法的jacobi矩阵
functionf=df1（x0）;
f=[2*x4*y
2*x-1];
以文件名df1.m保存。
步骤三：
编写主函数。
打开Editor编辑器输入以下语句：
%牛顿迭代法的主函数
x0=[0.80.7];
[x,n,data]=new_ton（x0）;
计算结果为'
迭代次数为'
n
以文件名new_main.m保存。
最速下降法原理
对于非线性方程组
令
如果给定一个初值
，我们希望找到一条路线每一次
迭代以后代价函数都会比原来小一些。
l称为步长因子
，
的不同，就构成了不同的下降算法。
如果取
就是所谓的最速下降法。
最速下降法是大范围收敛的h在某
出沿最速下降方向
下降的最快
symsxy
f1=x^2+2*y^2-2;
f2=x^2-y;
h=f1^2+f2^2
grad=[diff（h,x）,diff（h,y）];
grad=simple（grad）
以文件名tidu_fuhao.m保存。
编写梯度函数。
functionf=f1_tidu（x0）
f=[8*x^3+8*x*y^2-8*x-4*x*y
8*x^2*y+16*y^3-14*y-2*x^2]'
;
编写最速下降法的方法函数打开Editor编辑器，输入以下语句：
function[x,n,data]=zuisu（x0）
x1=x0-0.001*f1_tidu（x0）;
%迭代过程
tol）&
（n<
2000）
%0.001为步长因子
%data用来存放中间数据
以文件名zuisu.m保存。
步骤四：
%最速下降法的主函数
x0=[0.80.7];
[x,n,data]=zuisu（x0）;
用最速下降法的计算结果为'
由两种算法的运行结果截屏可看出，同样的精度要求（10-5）下，Newton法收敛速度远远好于最速下降法，我们知道最速下降法是线性收敛的。
上机作业4
1、用幂法与反幂法求矩阵A的按模最大、最小特征值与对应的特征向量。
幂法原理
当矩阵A满足一定条件时，在工程中可用幂法计算其主特征值（按模最大）及其特征向量。
矩阵A需要满足的条件为：
（1）
（2）存在n个线性无关的特征向量，设为
不全为0，则有
可见，当
越小时，收敛越快；
且当k充分大时，有
，对应的特征向量即是
。
设An有n个线性相关的特征向量v1，v2，…，vn，对应的特征值1，2，…，n，满足
|1|>
|2|…|n|
（3.2.1）
规范化
在实际计算中，若|1|>
1则|1ka1|，若|1|<
1则|1ka1|0都将停机。
须采用“规范化”的方法
，k=0，1，2，…
定理3.2-1任给初始向量
有，
function[t,y]=lpower（A,x0,eps,N）
A=[4,1,1,1;
1,3,-1,1;
1,-1,2,0;
1,1,0,2];
x0=[1;
1;
1];
eps=1e-5;
N=100;
k=1;
z=0;
y=x0./max（abs（x0））;
x=A*y;
b=max（x）;
ifabs（z-b）<
eps
t=max（x）;
return;
whileabs（z-b）>
eps&
&
k<
N
k=k+1;
z=b;
y=x./max（abs（x））;
[m,index]=max（abs（x））;
t=x（index）;
用幂法求矩阵的最大特征值及对应的特征向量'
最大特征值为：
t
最大特征值对应的特征向量为：
y
迭代次数：
k
反幂法理论
在工程计算中，可以利用反幂法计算矩阵按模最小特征值及其对应特征向量。
其基本理论如下，与幂法基本相同：
，可知，A和A-1的特征值互为倒数，求A按模最小特征值即求A-1的按模最大特征值,取倒数即为A的按模最小特征值所以算法基本相同，区别就是在计算
function[s,y]=invpower（A,x0,eps,n）
eps=10^（-5）;
N=100;
r=0;
[L,U]=lu（A）;
z=L\y;
x=U\z;
u=max（x）;
s=1/u;
whileabs（u-r）>
N
r=u;
s=1/x（index）;
用反幂法求矩阵的最小特征值及对应的特征向量'
最小特征值为：
s
最小特征值对应的特征向量为：
由上图可知，幂法与反幂法求解矩阵的特征值与特征向量是很方便的。
2、用Householder变换求矩阵A的QR分解，并用QR方法做3次迭代。
原理：
对任意一个非奇异矩阵（可逆矩阵）A，可以把它分解成一个正交阵Q和一个上三角阵的乘积，称为对矩阵A的QR分解，即A=QR。
如果规定R的对角元取正实数，这种分解是唯一的。
若A是奇异的，则A有零特征值。
任取一个不等于A的特征值的实数μ，则A‐μI是非奇异的。
只要求出A‐μI的特征值和特征向量就容易求出矩阵A的特征值和特征向量，所以假设A是非奇异的，不是一般性。
设A=A1，对A1作QR分解，得A1=Q1R1，交换该乘积的次序，得A2=R1Q1=，由于Q1正交矩阵，A1到A2的变换为正交相似变换，于是A1和A2就有相同的特征值。
一般的令A1=A，对k=1，2，3，…..
这样，可得到一个迭代序列{Ak}，这就是QR方法的基本过程。
QR算法的核心思想是：
对一个具有n个不相等的特征值的矩阵A，进行如下变换：
A=Q1R1;
A1=R1Q1
A1=Q2R2;
A2=R2Q2
…………
Ak=Qk+1Rk+1;
Ak+1=Rk+1Qk+1
其中Qi是正交矩阵，Ri是可逆的上三角矩阵。
由于：
Ai+1=QTi+1（Qi+1Ri+1）Qi+1=QTi+1AiQi+1
所以Ai+1和Ai相似，具有相同的特
升级会员 