Matlab非线性方程数值解法
实验目的
用Matlab实现非线性方程的二分法、不动点迭代法
实验要求
1. 给出二分法算法和不动点迭代算法
2. 用Matlab实现二分法
3. 用Matlab实现不动点迭代法
实验内容
（1）在区间[0,1]上用二分法和不动点迭代法求的根到小数点后六位。
（2）二分法的基本思想：逐步二分区间[a,b]，通过判断两端点函数值的符号，进一步缩小有限区间，将有根区间的长度缩小到充分小，从而，求得满足精度要求的根的近似值。
（3）不动点迭代法基本思想：已知一个近似根，构造一个递推关系（迭代格式），使用这个迭代格式反复校正根的近似值，计算出方程的一个根的近似值序列，使之逐步精确法，直到满足精度要求（该序列收敛于方程的根）。
实验步骤
（1）二分法算法与MATLAB程序（二分法的依据是根的存在性定理，更深地说是介值定理）。
MATLAB程序，
1 %二分法 2 %输入：f(x)=0的f(x)，[a,b]的a,b，精度ep 3 %输出：近似根root，迭代次数k 4 function [root,k]=bisect(fun,a,b,ep) 5 if nargin>3 6 elseif nargin<4 7 ep=1e-5;%默认精度 8 else 9 error(\'输入参数不足\');%输入参数必须包括f(x)和[a,b] 10 end 11 if fun(a)*fun(b)>0%输入的区间要求 12 root=[fun(a),fun(b)]; 13 k=0; 14 return; 15 end 16 k=1; 17 while abs(b-a)/2>ep%精度要求 18 mid=(a+b)/2;%中点 19 if fun(a)*fun(mid)<0 20 b=mid; 21 elseif fun(a)*fun(mid)>0 22 a=mid; 23 else 24 a=mid;b=mid; 25 end 26 k=k+1; 27 end 28 root=(a+b)/2; 29 end
二分法1
运行示例（并未对输出格式做控制，由于精度要求，事后有必要控制输出的精度）：
优化代码，减小迭代次数（在迭代前，先搜寻更适合的有根区间）
1 %二分法改良 2 %在一开始给定的区间中寻找更小的有根区间 3 %输入：f(x)=0的f(x)，[a,b]的a,b，精度ep 4 %输出：近似根root，迭代次数k 5 %得到的根是优化区间里的最大根 6 function [root,k]=bisect3(fun,a,b,ep) 7 if nargin>3 8 elseif nargin<4 9 ep=1e-5;%默认精度 10 else 11 error(\'输入参数不足\');%输入参数必须包括f(x)和[a,b] 12 end 13 %定义划分区间的分数 14 divQJ=1000; 15 %等分区间 16 tX=linspace(a,b,divQJ); 17 %计算函数值 18 tY=fun(tX); 19 %找到函数值的正负变化的位置 20 locM=find(tY<0); 21 locP=find(tY>0); 22 %定义新区间 23 if tY(1)<0 24 a=tX(locM(end)); 25 b=tX(locP(1)); 26 else 27 a=tX(locP(end)); 28 b=tX(locM(1)); 29 end 30 if fun(a)*fun(b)>0%输入的区间要求 31 root=[fun(a),fun(b)]; 32 k=0; 33 return; 34 end 35 k=1; 36 while abs(b-a)/2>ep%精度要求 37 mid=(a+b)/2;%中点 38 if fun(a)*fun(mid)<0 39 b=mid; 40 elseif fun(a)*fun(mid)>0 41 a=mid; 42 else 43 a=mid;b=mid; 44 end 45 k=k+1; 46 end 47 root=(a+b)/2; 48 end
二分法2
运行示例（同样没有控制输出）
明显地，迭代次数减小许多。
继续优化代码，以获得区间里所有的根（关键是对有根区间的判断与记忆）
1 %二分法改良 2 %在一开始给定的区间中寻找更小的有根区间 3 %输入：f(x)=0的f(x)，[a,b]的a,b，精度ep 4 %输出：近似根root，迭代次数k 5 %例子： 6 %fun=inline(\'x.^3-x.^2-4*x+4\');a=-3;b=3;ep=1e-5; 7 %fun=inline(\'x.^2-1\');a=-2;b=2;ep=1e-5; 8 %fun=inline(\'-x.^2+1\');a=-2;b=2;ep=1e-5; 9 %fun=inline(\'x.^2-3*x+2\');a=0;b=3;ep=1e-5; 10 function [root,k,x,y]=bisect4(fun,a,b,ep) 11 if nargin>3 12 elseif nargin<4 13 ep=1e-5;%默认精度 14 else 15 error(\'输入参数不足\');%输入参数必须包括f(x)和[a,b] 16 end 17 %定义划分区间的分数 18 divQJ=100; 19 %等分区间 20 tX=linspace(a,b,divQJ); 21 %计算函数值 22 tY=fun(tX); 23 %找到函数值的正负变化的位置 24 locM=find(tY<0); 25 locP=find(tY>0); 26 %得到的符号序列，找其中的间断点 27 z=0; 28 k=1; 29 for i=1:length(locM)-1%扫描-符号序列 30 %滚筒算法 31 if(i==1 && locM(1)~=1) 32 z=z+1; 33 x(z)=locM(i); 34 end 35 tmp=locM(i+1)-locM(i); 36 if(tmp~=1) 37 z=z+1; 38 x(z)=locM(i); 39 z=z+1; 40 x(z)=locM(i+1); 41 end 42 end 43 44 l=0; 45 for i=1:length(locP)-1%扫描+符号序列 46 %滚筒算法 47 if(i==1 && locP(1)~=1) 48 l=l+1; 49 y(l)=locP(i); 50 end 51 tmp=locP(i+1)-locP(i); 52 if(tmp~=1) 53 l=l+1; 54 y(l)=locP(i); 55 l=l+1; 56 y(l)=locP(i+1); 57 end 58 end 59 if(z==l-1) 60 z=z+1; 61 x(z)=locM(end); 62 elseif(z==l+1) 63 l=l+1; 64 y(l)=locP(i); 65 end 66 67 if(z==0 && l==0) 68 if(locM(end)<locP(1)) 69 x(1)=locM(end); 70 y(1)=locP(1); 71 else 72 x(1)=locM(1); 73 y(1)=locP(end); 74 end 75 z=1; 76 elseif(z==0) 77 x(1)=locM(1); 78 x(2)=locM(end); 79 z=2; 80 elseif(l==0) 81 y(1)=locP(1); 82 y(2)=locP(end); 83 z=2; 84 end 85 for i=1:z 86 a=tX(x(i)); 87 b=tX(y(i)); 88 if fun(a)*fun(b)>0%输入的区间要求 89 root=[fun(a),fun(b)]; 90 k=0; 91 return; 92 end 93 if(a>b) 94 tmp=b; 95 b=a; 96 a=tmp; 97 end 98 while abs(b-a)/2>ep%精度要求 99 mid=(a+b)/2;%中点 100 if fun(a)*fun(mid)<0 101 b=mid; 102 elseif fun(a)*fun(mid)>0 103 a=mid; 104 else 105 a=mid;b=mid; 106 end 107 k=k+1; 108 end 109 root(i)=(a+b)/2; 110 end 111 end
二分法3
运行示例（仍然没有控制输出，笔者忽然觉得自己中了月亮的毒！！！）
使用二分法解得在[0,1]上的根为0.090526。（终于在最后回答问题时做了<_> 控制输出的方法：在调用函数前，在控制台输入format long;得到一long型的数据，根据精度读取。）
（2）不动点迭代算法与MATLAB程序。
1 %不动点迭代法 2 %输入：fun函数,初始值x0，精度ep，最大迭代次数Nmax 3 %输出：近似根root，迭代次数k 4 function [root, k]=iterate(fun,x0,ep,Nmax) 5 if nargin<4 6 Nmax=500; 7 end 8 if nargin<3 9 ep=1e-5; 10 end 11 x=x0; 12 x0=x+2*ep; 13 k=0; 14 while abs(x0-x)>ep && k<Nmax 15 x0=x; 16 x=fun(x0); 17 k=k+1; 18 end 19 root=x; 20 if k==Nmax 21 warning(\'已达到迭代次数上限\'); 22 end 23 end
不动点迭代法
运行示例（对迭代方程的简单推导得到：）：
就没优化的二分法与不动点迭代法比较，明显地迭代法的迭代次数更少。
使用不动点迭代法，求得近似根0.090525。
小结
就上述两个程序而言，二分法的编程更具挑战。在优化二分法的程序过程中，发现二分法的关键就在于对所给定区间的讨论划分。得到一个能求解给定区间里的所有零点的程序的关键是（分类讨论的思想），对符号变化的记录与使用，换言之，解决这个问题的关键就在于分支结构的使用。
（不可否认，第三个二分法程序仍然有bug！！！待解决——有空再说。）