数值计算方法

实验一 非线性方程求根

一、实验目的

（1）熟悉二分法、牛顿迭代法求解非线性方程的计算过程；

（2）设计出相应的算法，编制相应的函数子程序。

二、实验内容

（1）画出二分法、牛顿迭代法的程序框图，编写并调试相应的函数子程序；

（2）分别用二分法、牛顿迭代法求解非线性方程在 (1, 2) 区间内的近似根，要求保留6位有效数字，并比较两种算法的收敛速度。

三、程序框图及实验步骤

3.1二分法程序框图：

3.2二分法实验步骤

3.2.1编写二分法程序

function [k,X]=bisection(f,a,b,Eps,N)

    

    f=inline('x^3+4*x^2-10');

    a=input('请输入解区间下限a:');

    b=input('请输入解区间上限b:');

    Eps=input('请输入精度eps:');

    N=input('请输入循环次数N:');

    for k =1:N

        A(k)=a;

        B(k)=b;

        ya=feval(f,a);

        yb=feval(f,b);

        temp=(a+b)/2;

        X(k)=temp;

        yt=feval(f,temp);

        F(k)=yt;

        if abs(yt)<Eps

            break;

        end

        if yt*ya<0

        a=a;

        b=temp;

        elseif yt*yb<0

            a=temp;

            b=b;

        end

    end

    sprintf(X);

end

3.2.2运行并输入初始值

3.2.3 得计算结果及运行时间

牛顿迭代程序框图

3.3牛顿迭代法实验步骤

3.3.1编写牛顿迭代法程序

function result= newdoniteration(x0,eps,N)

x0=input('请输入迭代初值x0:');

eps=input('请输入精度eps:');

N=input('请输入预计迭代次数:');

tic

y=inline('x^3 + 4*x^2 - 10');

dy=inline('3*x^2 + 8*x');

for k=1:N

    fx0=y(x0);

    dfx0=dy(x0);

    x1=x0-fx0/dfx0;

    if dfx0==0

        fprintf('导数为0');

        break    

    elseif abs(x1-x0)<eps

        fprintf('x*=%.5f\n',x0);

        fprintf('f(x*)=%.5f\n',f(x1));

        fprintf('迭代次数n=%d\n',k);            

        break

    else

        x0=x1;

    end    

end

toc

end

3.3.2运行并输入初始值

3..3.3 得计算结果及运行时间

3.4比较二分法与牛顿迭代法收敛速度

二分法迭代了18次，用时0.022190s

牛顿迭代法迭代了4次，用时0.005600s，由此可知牛顿迭代法得收敛速度高于二分法迭代速度。

实验二 线性方程组的迭代解法

1．实验目的

（1）熟悉雅可比迭代法、高斯-赛德尔迭代法求解线性方程组的计算过程；

（2）设计出相应的算法，编制相应的函数子程序。

2．实验内容

（1）画出雅可比迭代法、高斯-赛德尔迭代法的程序框图，编写并调试相应的函数子程序；

（2）分别用雅可比迭代法、高斯-赛德尔迭代法求解下述线性方程组，并比较两种算法的收敛速度。

3.程序框图及实验步骤

3.1雅可比迭代法

程序框图

3.1.1编写雅可比迭代法程序

function [x,n]=jacobi(A,b,x0,eps)

A=input('请输入系数矩阵A:');

b=input('请输入常数矩阵b:');

x0=input('请输入xi初始值矩阵x0:');

eps=input('请输入需要精度:');

D=diag(diag(A));%%求对角矩阵

L=-tril(A,-1);%%求下三角

U=-triu(A,1);%%求上三角

B=D\(L+U);

f=D\b;

x=B*x0+f;

n=1;

while norm(x-x0,'inf')>=eps %%用无穷范数去逼近

    Y=sprintf('% .4f',x),

    N=sprintf('%d',n),

    x0=x;

    x=B*x0+f;

    n=n+1;

End

3.1.2运行并初始化输入

3.1.3运行结果

3.2高斯——赛德尔方法

程序框图

3.2.1编写高斯赛德尔方法程序

function [x,n]=gauseidel(A,b,x0,eps)

A=input('请输入系数矩阵A:');

b=input('请输入常数矩阵b:');

x0=input('请输入xi初始值矩阵x0:');

eps=input('请输入需要精度:');

tic

D=diag(diag(A))%%求对角矩阵

L=-tril(A,-1);%%求下三角

U=-triu(A,1);%%求上三角

B=(D-L)\U;

f=(D-L)\b;

x=B*x0+f;

n=1;

while norm(x-x0,inf)>=eps %%用无穷范数去逼近

    Y=sprintf('% .4f',x),

    N=sprintf('%d',n),

    x0=x;

    x=B*x0+f;

    n=n+1;

end

toc

end

3.2.2运行程序及初始化输入

3.2.3运行结果

3.3比较雅可比方法与高斯赛德尔方法

由两方法得运算结果对比可知，在相同精度要求下，雅可比法迭代了18次，运算时间为0.011568s；高斯赛德尔迭代了9次，运算时间为0.008258s，因此，高斯赛德尔得收敛速度更快。

实验三 插值方法

1. 实验目的

（1）熟悉拉格朗日插值、牛顿插值等基本插值方法；

（2）设计出相应的算法，编制相应的函数子程序。

2．实验内容

（1）画出拉格朗日插值算法、牛顿插值算法的程序框图，编制并调试相应的函数子程序；

（2）已知函数四个点的数据如下表，试用拉格朗日插值确定函数在x=2.101，4.234处的函数值；

x	1.1	2.3	3.9	5.1
y	3.887	4.276	4.651	2.117

（3）已知试用牛顿插值公式求5^0.5的近似值。

3.程序框图与实验步骤

3.1拉格朗日程序框图

3.1.1编写拉格朗日程序

format long;

x0=[1.1,2.3,3.9,5.1];

y0=[3.887,4.276,4.651,2.117];

n=length(x0);

L=0;

l=1;

%x=str2double(input(''));

x=input('请输入需插值计算的x值:');

Tic

for i =1:n

    for j =1:n

        if i~=j

            l=l*(x-x0(j))/(x0(i)-x0(j));

            

        end

    end

    L=L+y0(i)*l;

    l=1;

end

fprintf('x=%.3f\n',L);

Toc

3.1.2运行及初始化输入

3.1.3运行结果

3.2牛顿插值方法

3.2.1牛顿插值程序框图

3.2.2编写牛顿插值程序

function [A,y]= newtoncz(X,Y,x)

x=input('请输入需计算的x:')

tic

X=[1,4,9];

Y=[1,2,3];

n=length(X); m=length(x);

for t=1:m

    z=x(t); A=zeros(n,n);A(:,1)=Y';

    s=0.0; y=0.0; c1=1.0;

    for  j=2:n

       for i=j:n

           A(i,j)=(A(i,j-1)- A(i-1,j-1))/(X(i)-X(i-j+1));

       end

    end

    for k=1:n

        p=1.0;

        for j=1:k-1

            p=p*(z-X(j));

        end

        s=s+A(k,k)*p;        

    end

    ss(t)=s;

end

    y=ss;

    A=[X',A];

    fprintf('y=%.4f\n',y);

toc

End

3.2.3运行程序并初始化输入

3.2.4运行结果

实验四 数值积分

1．实验目的

（1）熟悉复化梯形求积公式、复化辛浦生求积公式；

（2）设计出相应的算法，编制相应的函数子程序；

2．实验内容

（1）画出复化梯形求积算法、复化辛浦生求积算法的程序框图，编制并调试相应的函数子程序；

（2）分别用复化梯形公式、复化辛浦生公式计算定积分，取n=2，4，8，16。

3.程序框图及实验步骤

3.1复化梯形求积公式

3.1.1程序框图

3.1.2编写复化梯形求积公式程序

被积函数

function [fx] = Tfunction(x)

if x == 0

    fx = 1;

else

   fx = sin(x)/x;

end

end

复化梯形公式

function [result] = complexT( a, b,n)

a=input('请输入积分下限a:');

b=input('请输入积分上限b:');

n=input('请输入划分区间数n:');

tic

h = (b - a)/n;

result = 0;

for i = 1:n-1

    result = result + Tfunction(a+h*(i-1))+Tfunction(a+h*i);

end 

result = result * h / 2;

fprintf('积分结果=%.4f\n',result);

toc

end

3.1.3运行程序并初始化输入

3.1.4运行结果

n=2

n=4

n=8

n=16

3.2复化辛普森公式

3.2.1复化辛普森公式程序框图

3.2.2编写复化辛普森公式程序

被积函数

function [fx] = Tfunction(x)

if x == 0

    fx = 1;

else

   fx = sin(x)/x;

end

end

辛普森公式

function [result]=complexsimpson(a,b,n)

a=input('请输入积分下限a:');

b=input('请输入积分上限b:');

n=input('请输入划分区间数n:');

tic

h=(b-a)/n;

result=0;

xi=a:h:b;

for i =1:n

    result=result+4*Tfunction(xi(i)+0.5*h);

end

for k = 2:n

    result=result+2*Tfunction(xi(k));

end

result=(result+Tfunction(a)+Tfunction(b))*h/6;

fprintf('复化辛普森公式积分结果为：%.4f\n',result);

toc

End

3.2.3运行程序并初始化输入

    

3.2.4运行结果

n=2

3.3比较复化梯形公式与复化辛普森公式

由运行结果可知，复化梯形公式在划分区间数较少时，积分结果与实际结果的误差较大，收敛速度慢；复化辛普森公式在划分区间为2时，便比复化梯形公式划分16段的精度更高，收敛速度很快。