大学计算方法实验报告

数值实验报告

《计算方法》实验报告

实验题目

实验报告1：非线性方程组的求解···················P1~2

实验报告2：线性方程组解法·······················P3~4

实验报告3：Lagrange 插值多项式··················P5~7

姓名： 学号：

班级：

指导老师：

时间：

数值实验报告

1 专业

序号 日期

实验报告1：非线性方程组的求解

【实验目的】

1．用MATLAB来实践进行牛顿法的变形，即对牛顿法进行了修正，使其应用更为方便，掌握用MATLAB运用割线法求解非线性方程组。

2．运用MATLAB进行隐函数作图。

【实验内容】

[方法] 设a,b为迭代初值,求两点(a,f(a)) 与 (b,f(b)) 的连线(割线)与 x 轴的交点记为 c，再把迭代初值换成 b,c,重复计算.

[要求] 把下面程序复制为新的 M-文件,去掉开头的 %

再把 '?' 部分改写正确就是一个完整的程序,找前面一个例子试算

【解】在牛顿迭代公式中用差商代替导数。带入初值（a,f(a)）,（b,f(b)）,两点的连线与x轴的交点作为c，再把迭代初值换为b，c，重复计算。

【计算机求解】以y= x-exp(-x)为例初值a=0，b=1，误差不超过1.0*10^(-5)进行计算。

【程序如下】

function mysecant

f = inline('x-exp(-x)');

a = 0; b = 1;

delta = 10.^(-5); epsilon = 10.^(-5); max1 = 10;

[c,err,iter,yc] = secant(f,a,b,delta,epsilon,max1)

% ---------------------------------------------------------

function [c,err,iter,yc] = secant(f,a,b,delta,epsilon,max1)

% [c,err,iter,yc] = secant(f,a,b,delta,epsilon,max1)

% 输入: f 连续函数

% a,b 迭代初值

% delta,epsilon 容差

% max1 最大迭代次数

% 输出: c 近似根

% err 误差

% iter 迭代次数

% yc = f(c)

for k = 1:max1

ya = feval(f,a); % ya = f(a) 11kkkkxxxfxf)('kxf数值实验报告

2 yb = feval(f,b);

c = 1.1; % 割线与 x 轴交点的横坐标

err = abs(c-b); % 相邻两次迭代的误差

relerr = err/(abs(c)+eps); % 相对误差,eps 是matlab常数(机器精度)约为1e-16

% 为什么分母要加上一个小常数?

yc = feval(f,c);

if (err

break

end

a = b;

b = c;

end

iter = 8;

% -------------------------------------

【运行结果】

c = 1.1000

err = 0

iter = 8

yc = 0.7671

【结果分析】

迭代公式的收敛阶为要比简化牛顿迭代法公式收敛的快。

)()()()(111kkkkkkkxxxfxfxfxx618.1)51(21p数值实验报告

3

专业 序号 日期

实验报告2：线性方程组解法

【实验目的】

（1）高斯列主元消去法求解线性方程组的过程

（2）熟悉用迭代法求解线性方程组的过程

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

【实验内容】

Gauss列选主元消去法

首先完成下面Gauss列选主元的消去法程序,单独存为 gauss.m 文件(注意一定要与函数名相同）

然后找一个例子调用此程序验证是否正确(调用方法同Matlab内部函数调用完全一样)

【解】在消元的第k步，首先在第k列的下方所有元素中取绝对值最大的为主元，而后再进行消元。

【计算机求解】编写高斯消去法，运用matlab程序中矩阵来进行编写，计算结果。

【程序如下】

function x = gauss(A,b)

[n,n] = size(A);

x = zeros(n,1);

Aug = [A,b]; % 增广矩阵

for k = 1:n-1

[piv,r] = max(abs(Aug(k:n,k))); % 找列主元所在子矩阵的行r

r = r + k - 1; % 列主元所在大矩阵的行

if r>k

temp=Aug(k,:);Aug(k,:)=Aug(r,:);Aug(r,:)=temp;

% %对Aug实施行交换

end

if Aug(k,k)==0, error('对角元出现0'), end % 程序遇到error会中断执行并显示其中的提示内容

%% 把增广矩阵消元成为上三角 数值实验报告

4 for p = k+1:n

mult = Aug(p,k)/Aug(k,k); % 消元乘子

Aug(p,k:n+1) = Aug(p,k:n+1)-Aug(k,k:n+1)*mult;

end

end

% % 解上三角方程组

A = Aug(:,1:n); b = Aug(:,n+1);

x(n) = b(n)/A(n,n);

for k = n-1:-1:1

x(k) = b(k);

for p=n:-1:k+1

x(k)=x(k)-A(k,p)*x(p);

end

x(k)=x(k)/A(k,k);

end

【运行结果】

A=[4 2 1;2 -5 2;1 2 6];b=[7;-1;9];

gauss(A,b)

ans =

1

1

1

【结果分析】

这个题目如果使用MATLAB自带的求解函数fzero来求得更精确的解，所得到的结果为0.7391。观察所得到的结果并进行比较，该方法的求解效果是较好的，达到了所要求的精度要求。说明此方法的效果很好。

数值实验报告

5

专业 序号 日期

实验报告3：Lagrange 插值多项式

【实验目的】

通过matlab计算作图，体会节点的增多虽然会使插值多项式P（x）在更多的地方与f（x）相等，但是两个插值之间P（x）不一定能很好的逼近f(x)，有时甚至差别很大，即会产生龙格现象。

【实验内容】

Runge现象的演示

在区间[-1,1]上给定函数f(x)=1./(1+25*x.^2 设将区间[-1,1]分为n等分，作出f(x)在此区间产生的龙格现象。

【实验过程】

function y = mylagpoly(X,Y,x)

f = inline('1./(1+25*x.^2)');

n = 11;

X = linspace(-1,1,n);

Y = f(X);

x = -1 : 0.01 : 1;

n = length(X);

y = zeros(size(x));

for i = 1:n

t=1;

for j = 1:n

if j ~= i

t = t.*((x-X(j))/(X(i)-X(j)));

end

end

y = y+t.*Y(i);

end

plot(x,f(x),'r',X,Y,'p',x,y,'k')

title('Runge现象')

legend('y=1/(1+25*x^2)''插值点 ','等分的10次插值多项式',0)

function [C,L] = interp_lag(X,Y)

n = length(X);

L = zeros(n,n);

for j = 1:n