实验报告一、求方程f(x)=x^3-sinx-12x+1的全部根, ε=1e -6 1、 用一般迭代法; 2、 用牛顿迭代法;并比较两种迭代的收敛速度。
一、首先,由题可求得:12cos 3)(2'--=x x x f . 其次,分析得到其根所在的区间。
① 令()0=x f ,可得到x x x sin 1123=+-.② 用一阶导数分析得到1123+-x x 和x sin 两个函数的增减区间;再用二阶导数分析得到两个函数的拐点以及凹凸区间.③ 在直角坐标轴上描摹出01123=+-x x 和0sin =x 的图,在图上可以看到他们的交点,然后估计交点所在的区间,即是所要求的根的区间。
经过估计,得到根所在的区间为[]3,4--,[]1,0和[]4,3.1、 一般迭代法 (1)算法步骤:设ε为给定的允许精度,迭代法的计算步骤为:① 选定初值0x .由()0=x f 确定函数()x g ,得等价形式()x g x =. ② 计算()0x g .由迭代公式得()01x g x =.③ 如果ε≤-01x x ,则迭代结束,取1x 为解的近似值;否则,用1x 代替0x ,重复步骤②和步骤③. (2)程序代码: ① 在区间[]3,4--内, 代码: clcx0=-3.5; %初值0xiter_max=100; %迭代的最大次数 ep=1e-6; %允许精度 ε k=0;while k<=iter_max %k 从0开始到iter_max 循环 x1=(sin(x0)+12*x0-1).^(1/3); %代入0x ,算出1x 的值 if abs(x1-x0)<ep %01x x -与允许精度作比较 break ; %条件ε≤-01x x 成立,跳出循环 endx0=x1; %条件ε≤-01x x 不成立,用1x 代替0x k=k+1; %k 加1 endx_star=x1, iter=k %1x 为解的近似值,iter 为迭代次数 运行结果:x_star = -3.4101 ;iter =14②在区间[]1,0内,代码: clcx0=0.5; %初值0xiter_max=100; %迭代的最大次数 ep=1e-6; %允许精度ε k=0;while k<=iter_max %k 从0开始到iter_max 循环 x1=(1/12)*(x0.^3-sin(x0)+1); %代入0x ,算出1x 的值if abs(x1-x0)<ep %01x x -与允许精度作比较 break ; %条件ε≤-01x x 成立,跳出循环 endx0=x1; %条件ε≤-01x x 不成立,用1x 代替0x k=k+1; %k 加1 endx_star=x1, iter=k %1x 为解的近似值,iter 为迭代次数 结果:x_star = 0.07696;iter =6③在区间[]4,3内, 代码: clcx0=3.5; %初值0xiter_max=100; %迭代的最大次数 ep=1e-6; %允许精度ε k=0;while k<=iter_max %k 从0开始到iter_max 循环 x1=(sin(x0)+12*x0-1).^(1/3); %代入0x ,算出1x 的值 if abs(x1-x0)<ep %01x x -与允许精度作比较 break ; %条件ε≤-01x x 成立,跳出循环 endx0=x1; %条件ε≤-01x x 不成立,用1x 代替0x k=k+1; %k 加1 endx_star=x1, iter=k %1x 为解的近似值,iter 为迭代次数 运行结果:x_star = 3.4101 ;iter =102、 牛顿迭代法 (1)算法步骤:① 选定初值0x ,计算()0x f ,()0'x f .② 按公式()()k k k k x f x f x x '1-=+迭代,得新的近似值1+k x ,并计算()1+k x f ,()1'+k x f . ③ 对于给定的允许精度ε,如果ε≤-+k k x x 1,则终止迭代,取1*+≈k x x ;否则,1+=k k ,在转回步骤②计算.(2)程序代码: ①在区间[]3,4--内, clcx1=-3.5; %初值1x k=0;while k<=100 %k 从0开始到100循环 x0=x1; %将初值1x 赋给0xf0=x0.^3-sin(x0)-12*x0+1; %计算()0x f f1=3*x0.^2-cos(x0)-12; %计算()0'x f x1=x0-f0/f1; %计算得到新的近似值1xif abs(x1-x0)< 1.0e-6 %01x x -与允许精度作比较 break ; %条件ε≤-01x x 成立,跳出循环 endk=k+1; %条件ε≤-01x x 不成立,k 加1endx_star=x1, iter=k %1x 为解的近似值,iter 为迭代次数 运行结果:x_star = -3.4911;iter =2②在区间[]1,0内, clcx1=0.5; %初值1x k=0;while k<=100 %k 从0开始到100循环 x0=x1; %将初值1x 赋给0xf0=x0.^3-sin(x0)-12*x0+1; %计算()0x f f1=3*x0.^2-cos(x0)-12; %计算()0'x f x1=x0-f0/f1; %计算得到新的近似值1xif abs(x1-x0)< 1.0e-6 %01x x -与允许精度作比较 break ; %条件ε≤-01x x 成立,跳出循环 endk=k+1; %条件ε≤-01x x 不成立,k 加1 endx_star=x1, iter=k %1x 为解的近似值,iter 为迭代次数 运行结果:x_star =0.07696 ;iter =3③在区间[]4,3内, clcx1=3.5; %初值1x k=0;while k<=100 %k 从0开始到100循环 x0=x1; %将初值1x 赋给0xf0=x0.^3-sin(x0)-12*x0+1; %计算()0x f f1=3*x0.^2-cos(x0)-12; %计算()0'x f x1=x0-f0/f1; %计算得到新的近似值1xif abs(x1-x0)< 1.0e-6 %01x x -与允许精度作比较 break ; %条件ε≤-01x x 成立,跳出循环 endk=k+1; %条件ε≤-01x x 不成立,k 加1 endx_star=x1, iter=k %1x 为解的近似值,iter 为迭代次数 运行结果:x_star =3.4101;iter =33、运行结果:4、结果分析:从这题的结果可以看出,牛顿迭代法的迭代速度比一般迭代法的速度要快,牛顿法的迭代次数比较少。
例如,求在[]3,4--的根时,一般迭代法迭代了14次(iter =14),而牛顿法只用了2次(iter =2)。
总而言之,一般迭代法是一种逐次逼近的方法,具有原理简单、编制程序方便等优点,但存在是否收敛和收敛速度的快慢等问题。
牛顿迭代法是一种特殊的迭代法,用于求方程的单根时具有2阶收敛。
因此是一种求非线性方程解的好方法,还可以用于求重根和复根,而且可以推广到求非线性方程组的解.二、解方程组直接法1、已知42158721048361261120A ⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥⎣⎦ 对矩阵A 做LU 分解。
2、用追赶法解下述方程组Ax b =,并给出n=10的结果,其中⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡=2112112112 A ,⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡----=5557 b1、杜利特尔分解法(直接三角分解法):设方程组b Ax =的系数矩阵A 的各阶主子式()()n k A k ,,2,10det =≠,则可知该方程组存在唯一的杜利特尔分解:LU A =,其中⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=-11112121nn n n l l l l L ,⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=nn n n u u u u u u U 22211211. 我们记矩阵A 的第i 行第j 列元素为ij a ,则矩阵U 的第一行和L 的第一列可由公式(1.1)求出:⎪⎩⎪⎨⎧====.,,3,2,,,,2,1,111111n i ua l n j a u i i j j (2.1) 再由公式(2.2)求出U 的第k 行()n k ,,2 =、L 的第k 列()1,,2-=n k 元素:⎪⎪⎪⎩⎪⎪⎪⎨⎧++=-=+=-=∑∑-=-=.,,2,1,,,,1,,1111n k k i u u l a l n k k j u l a u kk k p pkip ik ik k p pj kp kj kj (2.2)(1)算法步骤:由于系数矩阵A 的各阶主子式()()n k A k ,,2,10det =≠,则可知该方程组存在唯一的杜利特尔分解:LU A =。
① 利用公式(2.1)和(2.2)求出U 和L 的元素ki u 和ik l 。
② 求解单位下三角形方程组b Ly =,即按公式⎪⎩⎪⎨⎧=-==∑-=1111,,3,2,,k j jkj k k n k y l b y b y (2.3)求出n y y y ,,,21 。
④ 求解上三角形方程组y Ux =,即按公式⎪⎪⎪⎩⎪⎪⎪⎨⎧--=-==∑+=1,,2,1,,1n n k u x u y x u y x kk nk j jkj k k nnn n (2.4) 可求得方程组的解:11,,,x x x n n -。
(2)程序代码:function [L,U]=LU(A)A[n,n]=size(A); %定义A 为n n ⨯阶矩阵 L=zeros(n,n); %矩阵L 初始化 U=zeros(n,n); %矩阵U 初始化 for i=1:nL(i,i)=1; %矩阵L 的对角线元素为1endfor k=1:n for j=k:nU(k,j)=A(k,j)-sum(L(k,1:k-1).*U(1:k-1,j)');%求出U 元素ki uendfor i=k+1:nL(i,k)=(A(i,k)-sum(L(i,1:k-1).*U(1:k-1,k)'))/U(k,k);%求出L 的元素ik lend endA=[4 2 1 5;8 7 2 10;4 8 3 6;12 6 11 20]; %矩阵A [L,U]=LU(A) %LU A =2、追赶法设所求方程组为⎪⎪⎪⎪⎩⎪⎪⎪⎪⎨⎧=+=++=++=++=+-------n n n n n n n n n n n n r x b x a r x c x b x a r x c x b x a r x c x b x a r x c x b 11111213433323232221212111 , (2.5)(1)算法步骤:① 由方程组的第1个方程,将未知元1x 用2x 表示为12111b x c r x -=记111b r u =,111b cv =,得 2111x v u x -=。