如何用二分法求解非线性方程

 

matlab是做数值运算尤其是矩阵运算非常牛逼的工具。并且可以嵌入到其他的各种开发工具当中，功能非常的强大！

毫不夸张的说，除了生孩子。matlab什么都能干。哈哈

好废话不多说，直接教大家如何用二分法求解非线性方程。

 

实际上，一般情况下，我们可以用二分法，迭代法，牛顿法，弦截法。这几种方法各有优缺点，在这里简单的分析一下：

二分法

其实质是对连续函数采用不断地减半缩小区间的做法来达到逼近真实值，从而求解出该方程的解。

 仅要求f（x）在[a，b]上连续,且f（a）*f（b）<0

优点是：算法比较简单

缺点是：收敛速度比较慢

 

迭代法：

所谓迭代法就是采用就是用某种收敛于所给问题的精确解的极限过程来逐步逼近的一种计算方法，从而可以用有限的步骤算出精确解的具有指定精度近似解。

对于一般形式的f（x）=0，先化简成x=g（x），再从每一个数列x0出发，做序列Χ_k+1=f(Χ_k),k=0,1,2....当k取极限的时候，得出解。

优点：收敛速度比二分法快，程序算法的编写难度也还好

缺点：迭代函数不唯一，可能收敛，可能发散。一般要选择能使g^`(x)<1的。

 

牛顿法：

牛顿法与迭代法有点相似，但是它是以做切线的方式，即求导来不断逼近真实值的。

优点：收敛速度很快

缺点：程序算法的编写比较有难度，并且由于需要求导，对于一下比较复杂的非线性方程。如超越方程等...不好求导，计算量大！

 

弦截法：

由于我们知道牛顿法的收敛速度非常快，但是为了解决那些难于直接求导的非线性方程，我们使用弦截法。他是利用说f^`(x)=(y2-y1)/(x2-x1),从而到达免去复杂的求导运算。

优点：收敛速度比较快，实质是牛顿法的一种变形

缺点：算法复杂程度还能接受

 

相信现在大家对这集中基本处理非线性方程的算法有了基础的了解，建议去看看数值计算的书。这样子能理解得更深刻。

正是回来主题，用二分法来解决光纤传播过程中传播常数的问题。

实际上我对光纤传导问题的理解也是一知半解。只知道说，中间介质折射率大，两边介质折射率小的时候。光以一定的角度入射（布如斯特角），就能发生却放射，一直传播一下。当然啦信号肯定还是会有衰减的。

算法的流程图设计如下：

然后有两个m文件，一个是描述方程的，一个是具体算法的实现过程。

方程的m文件：

 %An M-file which defines a function needed be solved
function f=fun(beita,N)
%变量的赋值
nameda=1.55*10^(-6);
k=2*pi/nameda;
n1=1.5;
n2=1.45;
n3=1.45;
A=2*10^(-6);
%中间代换量
p=sqrt((beita)^2-(n2*k)^2);
q=sqrt((n1*k)^2-(beita)^2);
r=sqrt((beita)^2-(n3*k)^2);
%方程
f=4*A*q-2/tan(p/q)-2/tan(r/q)-2*N*pi;

算法的m文件：

%An M-file to calculate the function by dichotomy
%变量定义
disp('利用二分法求传播常数beita');
n1=1.55;
n2=1.45;
nameda=1.55*10^(-6);
k=2*pi/nameda;
b=n1*k;
a=n2*k;
h=(b-a)/100;
%二分法具体实现过程
eps=input('请输入eps的值：');
    a1=a;
    b1=a1+h;
    N=0;
if N<10
    while b1<b
        while fun(a1,N)*fun(b1,N)>0
            a1=b1;
            b1=a1+h;
        end

        while abs(b1-a1)>eps
            x=0.5*(a1+b1);
            f2=fun(x,N);

            if f2==0 
                continue, end   
            
            if fun(a1,N)*fun(x,N)<0
                b1=x;
            else
                a1=x;
            end
        end  

        fid = fopen('data.txt', 'at+');
        fprintf(fid, 'N=%d\t\tbeita=%d \n',N,x);
        fclose(fid);
        a1=x+h/10;
        b1=a1+h;
        N=N+1;
    end
end
disp('计算完毕，请查阅data.txt文件');

到这里，利用二分法求解光纤传导问题的求解过程就结束了。

下次，我弄个用弦截法求解的上来。