数值分析上机作业

合集下载

数值分析第二次上机

题目：在区间[-1,1]上分别取n=10、20用两组等距节点对龙格函数22511)(x x f +=做多项式插值及三次样条插值，对每个n 值，分别画出插值函数及)(x f 的图形程序：

function y1=lagrange(x0,y0,x1)

n=length(x0);

syms x;

for k=1:n

l(k)=x/x;

for p=1:n

if p~=k

l(k)=l(k)*(x-x0(p))/(x0(k)-x0(p));

end

z=0;

for k=1:n

z=z+l(k)*y0(k);

end

y1=subs(z,x,x1);

程序：n=10时

n=10;

x0=-1:2/n:1;

y0=1./(1+25*x0.^2);

x=-1:.001:1;

y1=lagrange(x0,y0,x);

y=1./(1+25*x.^2);

plot(x,y,x,y1,'-.',x0,y0,'p');

legend('Runge Function','插值函数','插值节点');

title('n=10时的Lagrange 插值的龙格现象');

xlabel('x');

ylabel('y');

结果：

程序：n=20时

n=20;

x0=-1:2/n:1;

y0=1./(1+25*x0.^2);

x=-1:.001:1;

y1=lagrange(x0,y0,x);

y=1./(1+25*x.^2);

plot(x,y,x,y1,'-.',x0,y0,'p');

legend('Runge Function','插值函数','插值节点'); title('n=20时的Lagrange插值的龙格现象'); xlabel('x');

数值分析第二次上机作业实验报告

一．实验任务

用MA TLAB 语言编写连续函数最佳平方逼近的算法程序（函数式M 文件）。并用此程序进行数值试验，写出实验报告。

二．实验方法

最佳平方逼近方法采用基于正交多项式的最佳平方逼近，选择Lengendre 多项式做基。计算组合系数时，函数的积分采用变步长复化梯形求积法。

三．程序功能和使用说明

1.采用基于正交多项式的最佳平方逼近，选择Lengendre 多项式做基

利用递推关系

0112()1,()()(21)()(1)()/2,3,.....

n n n P x P x x

P x n xP x n P x n n --===---⎡⎤⎣⎦=

可构造出用户需要的任意次数的最佳平方逼近多项式。

2. 用M 文件建立数学函数，实现程序通过修改建立数学函数的M 文件以适用不同的被逼近函数。

3.已经考虑一般的情况]1,1[],[)(+-≠∈b a x f ，程序有变量代换的功能。

4.计算组合系数时，函数的积分采用变步长复化梯形求积法

5.可根据需要，求出二次、三次、。。。最佳平方逼近函数）x s (。

6.最后作出逼近函数）x s (和被逼近函数)(x f 的曲线图可进行比较，分别用绘图函数plot 和fplot 绘图。

7.在matlab 的命令窗口，输入[c,sx]=leastp(@func1,a,b,n)，func1是被逼近函数，b 和a 分别是逼近函数的上、下区间，n 为最佳平方逼近的次数，可为任意次数。

四．程序代码（含注释）

1. 最佳平方逼近主函数

function [c,sx]=leastp(func,a,b,n)

东南大学数值分析上机作业汇总

东南大学数值分析上机作业

汇总

-标准化文件发布号：（9456-EUATWK-MWUB-WUNN-INNUL-DDQTY-KII

数值分析上机报告

院系：

学号：

姓名：

作业1、舍入误差与有效数 (1)

1、函数文件cxdd.m (1)

2、函数文件cddx.m (1)

3、两种方法有效位数对比 (1)

4、心得 (2)

作业2、Newton迭代法 (2)

1、通用程序函数文件 (3)

2、局部收敛性 (4)

（1）最大δ值文件 (4)

（2）验证局部收敛性 (4)

3、心得 (6)

作业3、列主元素Gauss消去法 (7)

1、列主元Gauss消去法的通用程序 (7)

2、解题中线性方程组 (7)

3、心得 (9)

作业4、三次样条插值函数 (10)

1、第一型三次样条插值函数通用程序： (10)

2、数据输入及计算结果 (12)

作业1、舍入误差与有效数

设∑

=-=N

j N j S 2

，其精确值为⎪⎭

⎫ ⎝⎛---1112321N N . （1）编制按从小到大的顺序1

131121222-⋅

⋅⋅+-+-=N S N ，计算N S 的通用程序；

（2）编制按从大到小的顺序()1

11111222-⋅⋅⋅+--+-=N N S N ，计算N S 的通用程序；

（3）按两种顺序分别计算642101010,,S S S ,并指出有效位数；（4）通过本上机你明白了什么？程序：

1、函数文件cxdd.m

function S=cxdd(N) S=0; i=2.0;

while (i<=N) S=S+1.0/(i*i-1); i=i+1;

数值分析上机作业(总)

数值分析上机实验

一、解线性方程组直接法（教材49页14题）

追赶法程序如下：

function x=followup(A,b)

n = rank(A);

for(i=1:n)

if(A(i,i)==0)

disp('Error: 对角有元素为0');

return;

end

end;

d = ones(n,1);

a = ones(n-1,1);

c = ones(n-1);

for(i=1:n-1)

a(i,1)=A(i+1,i);

c(i,1)=A(i,i+1);

d(i,1)=A(i,i);

end

d(n,1) = A(n,n);

for(i=2:n)

d(i,1)=d(i,1) - (a(i-1,1)/d(i-1,1))*c(i-1,1);

b(i,1)=b(i,1) - (a(i-1,1)/d(i-1,1))*b(i-1,1);

end

x(n,1) = b(n,1)/d(n,1);

for(i=(n-1):-1:1)

x(i,1) = (b(i,1)-c(i,1)*x(i+1,1))/d(i,1);

end

主程序如下：

function zhunganfa

A=[2 -2 0 0 0 0 0 0;-2 5 -2 0 0 0 0 0;0 -2 5 -2 0 0 0 0;0 0 -2 5 -2 0 0 0;0 0 0 -2 5 -2 0 0;0 0 0 0 -2 5 -2 0;0 0 0 0 0 -2 5 -2;0 0 0 0 0 0 -2 5];

b=[220/27;0;0;0;0;0;0;0];

x=followup(A,b)

数值分析大作业

数值分析上机作业（一）

一、算法的设计方案

1、幂法求解λ1、λ501

幂法主要用于计算矩阵的按模最大的特征值和相应的特征向量，即对于

|λ1|≥|λ2|≥.....≥|λn|可以采用幂法直接求出λ1，但在本题中λ1≤λ2≤……≤λ501，我们无法判断按模最大的特征值。但是由矩阵A的特征值条件可知|λ1|和|λ501|之间必然有一个是最大的，通过对矩阵A使用幂法迭代一定次数后得到满足精度ε=10−12的特征值λ0，然后在对矩阵A做如下的平移：

B=A-λ0I

由线性代数(A-PI)x=(λ-p)x可得矩阵B的特征值为：λ1-λ0、λ2-λ0…….λ501-λ0。对B矩阵采用幂法求出B矩阵按模最大的特征值为λ∗=λ501-λ0，所以λ501=λ∗+λ0，比较λ0与λ501的大小，若λ0>λ501则λ1=λ501，λ501=λ0；若λ0

求矩阵M按模最大的特征值λ的具体算法如下：

任取非零向量u0∈R n

ηk−1=u T(k−1)∗u k−1

y k−1=u k−1η

k−1

u k=Ay k−1

βk=y T

k−1u k

(k=1,2,3……)

当|βk−βk−1|

|βk|

≤ε=10−12时，迭终终止，并且令λ1=βk

2、反幂法计算λs和λik

由已知条件可知λs是矩阵A 按模最小的特征值，可以应用反幂法直接求解出λs。

使用带偏移量的反幂法求解λik，其中偏移量为μk=λ1+kλ501−λ1

(k=1,2,3…39)，构造矩阵C=A-μk I，矩阵C的特征值为λik−μk，对矩阵C使用反幂法求得按模

最小特征值λ0，则有λik=1

(完整版)哈工大-数值分析上机实验报告

实验报告一

题目：非线性方程求解

摘要：非线性方程的解析解通常很难给出，因此线性方程的数值解法就尤为重要。本实验采用两种常见的求解方法二分法和Newton法及改进的Newton法。

前言：（目的和意义）

掌握二分法与Newton法的基本原理和应用。

数学原理：

对于一个非线性方程的数值解法很多。在此介绍两种最常见的方法：二分法和Newton法。

对于二分法，其数学实质就是说对于给定的待求解的方程f(x)，其在[a,b]上连续，f(a)f(b)<0，且f(x)在[a,b]内仅有一个实根x*，取区间中点c，若，则c恰为其根，否则根据f(a)f(c)<0是否成立判断根在区间[a,c]和[c,b]中的哪一个，从而得出新区间，仍称为[a,b]。重复运行计算，直至满足精度为止。这就是二分法的计算思想。

Newton法通常预先要给出一个猜测初值x0，然后根据其迭代公式

产生逼近解x*的迭代数列{x k}，这就是Newton法的思想。当x0接近x*时收敛很快，但是当x0选择不好时，可能会发散，因此初值的选取很重要。另外，若将该迭代公式改进为

其中r为要求的方程的根的重数，这就是改进的Newton法，当求解已知重数的方程的根时，在同种条件下其收敛速度要比Newton法快的多。

程序设计：

本实验采用Matlab的M文件编写。其中待求解的方程写成function的方式，如下

function y=f(x);

y=-x*x-sin(x);

写成如上形式即可，下面给出主程序。

二分法源程序：

clear

%%%给定求解区间

b=1.5;

a=0;

数值分析上机实验报告

数

值

分

析

上

机

实

验

理学院

11级统计01班

41108030125

鲁庆

实验报告一

一．实验名称

误差与误差估计

二．实验目的

掌握数值运算的误差估计方法

三．数学原理 1.绝对误差(*)e x

设某一量的准确值为x ，近似值为x*，则x*与x 之差叫做近似值x*的绝对误差(简称误差)，记为*(*)*e e x x x ==- 2.绝对误差限

适当小的正数，使|(*)||*|*e x x x ε=-≤则称*ε为近似值 x * 的绝对误差限。(有时用

*x x ε*

=±表示近似值x *的精度或准确值的所在范围。

3.相对误差(*)r e x

绝对误差与准确值之比*

(*)*(*),0r r e x x x

e e x x x x

-===≠称为x *的相对误差。

4.相对误差限(*)r x ε

若指定一个适当小的正数 (*)r x ε，使|(*)|

|(*)|(*)||

r r e x e x x x ε=≤则称(*)r x ε为近似值 x *的相对误差限。

5.有效数字

若近似值x*的绝对误差限是某一位的半个单位，该位到x*的第一位非零数字一共有n 位，则称近似值x*有n 位有效数字，或说x*精确到该位。 6.绝对误差的运算：

)()()(2121x x x x εεε+=± )()()(122121x x x x x x εεε+≈

122121

+=x x x x x x x )()()(

εεε (f(x))()(x)f x εε'≈

四．实验内容

1. 计算

I n

=e 1

-⎰1

e x 2

dx (n=0,1,...)并估计误差。

数值分析上机实验

1 绪论 (1)

2 实验题目(一) (2)

2.1 题目要求 (2)

2.2 NEWTON插值多项式 (3)

2.3 数据分析 (4)

2.3.1 NEWTON插值多项式数据分析 (4)

2.3.2 NEWTON插值多项式数据分析 (6)

2.4 问答题 (6)

2.5 总结 (7)

3 实验题目(二) (8)

3.1 题目要求 (8)

3.2 高斯-塞德尔迭代法 (8)

3.3 高斯-塞德尔改进法—松弛法 (9)

3.4 松弛法的程序设计与分析 (9)

3.4.1 算法实现 (9)

3.4.2 运算结果 (9)

3.4.3 数据分析 (11)

4 实验题目(三) (13)

4.1 题目要求 (13)

4.2 RUNGE-KUTTA 4阶算法 (13)

4.3 RUNGE-KUTTA 4阶算法运算结果及数值分析 (14)

总结 (16)

附录A (17)

1绪论

数值分析是计算数学的一个主要部分，它主要研究各类数学问题的数值解法，以及分析所用数值解法在理论上的合理性。实际工程中的数学问题非常复杂，所以往往需要借助计算机进行计算。运用数值分析解决问题的过程：分析实际问题，构建数学模型，运用数值计算方法，进行程序设计，最后上机计算求出结果。

数值分析这门学科具有面向计算机、可靠的理论分析、好的计算复杂性、数值实验、对算法进行误差分析等特点。

本学期开设了数值分析课程，该课程讲授了数值分析绪论、非线性方程的求解、线性方程组的直接接法、线性方程组的迭代法、插值法、函数逼近与曲线拟合、数值积分和数值微分、常微分方程初值问题的数值解法等内容。其为我们解决实际数学问题提供了理论基础，同时我们也发现课程中很多问题的求解必须借助计算机运算，人工计算量太大甚至无法操作。所以学好数值分析的关键是要加强上机操作，即利用计算机程序语言实现数值分析的算法。本报告就是基于此目的完成的。

硕士研究生数值分析上机实验

数值分析上机实验

1、已知A 与b

12.38412 2.115237 -1.061074 1.112336 -0.1135840.718719 1.742382 3.067813 -2.031743 2.11523719.141823

-3.125432 -1.012345 2.189736

1.563849

-0.784165 1.112348 3.123124 -1.061074 -3.125A =43215.567914 3.123848 2.031454 1.836742-1.056781 0.336993 -1.010103 1.112336 -1.012345 3.12384827.108437 4.101011-3.741856 2.101023 -0.71828 -0.037585 -0.113584

2.189736 2.031454 4.10101119.8979180.431637-

3.111223 2.121314 1.784317 0.718719 1.563849 1.836742 -3.741856 0.4316379.789365-0.103458 -1.103456 0.238417 1.742382 -0.784165 -1.056781 2.101023-3.111223-0.1034581

4.7138465 3.123789 -2.213474 3.067813 1.112348 0.336993-0.71828 2.121314-1.103456 3.12378930.719334 4.446782 -2.031743 3.123124 -1.010103-0.037585 1.7843170.238417-2.213474 4.44678240.00001[ 2.1874369 33.992318 -2

数值分析上机第四次作业

实验项目：共轭梯度法求解对称正定的线性方程组

实验内容：用共轭梯度法求解下面方程组

(1) 123421003131020141100155x x x x -⎛⎫⎛⎫⎛⎫ ⎪ ⎪ ⎪--- ⎪ ⎪ ⎪= ⎪ ⎪ ⎪-- ⎪ ⎪ ⎪-⎝⎭⎝⎭

⎝⎭ 迭代20次或满足()(1)

1110k k x x --∞-<时停止计算。

(2) Ax b =,A 是1000阶的Hilbert 矩阵或如下的三对角矩阵，

A[i,i]=4，A[i,i-1]=A[i-1,i]=-1，i=2,3,..,n

b[1]=3, b[n]=3, b[i]=2，i=2,3,…,n-1

迭代10000次或满足()()710k k r b Ax -=-≤时停止计算。

（3）*考虑模型问题，方程为

222222(),(,)(0,1)(0,1)(0,)1,(1,),01

(,0)1,(,1),01xy y x u u x y e x y D x y

u y u y e y u x u x e x ∂∂+=+∈=⨯∂∂==≤≤==≤≤

用正方形网格离散化，若取1/,10h N N ==，得到100n =的线性方程组，并用共轭梯度法（CG 法）求解，并对解作图。

实验要求：迭代初值可以取01(,1,...,)ij u i j N ==，计算到32||||10k r -≤停止．本

题有精确解(,)xy u x y e =，这里k u 表示以k ij u 为分量的向量，u 表示在相应点(,)i j 上

取值作为分量的向量．

实验一：

（1）

高等数值分析上机作业

,3 ,3 ,3
) ) )
a1 a2 a3
( ( (
f f f
,1 ,2 ,3
) ) )
其中0 (x) 1,1(x) x,2 (x) x2,3(x) x3, f (x) arcsin x.
1
程序代码： a=zeros(4); syms x a(1,1)=int('1',-1,1);a(1,2)=int('x',-1,1); a(1,3)=int('x^2',-1,1);a(1,4)=int('x^3',-1,1); a(2,2)=int('x^2',-1,1);a(2,3)=int('x^3',-1,1); a(2,4)=int('x^4',-1,1);a(3,3)=int('x^4',-1,1); a(3,4)=int('x^5',-1,1);a(4,4)=int('x^6',-1,1); a(2,1)=a(1,2);a(3,1)=a(1,3);a(4,1)=a(1,4); a(3,2)=a(2,3);a(4,2)=a(2,4);a(4,3)=a(3,4); a B=zeros(4,1); h=0.01;x=-1:h:1;y=asin(x).*x.^3; B(1,1)=quad('asin(x)',-1,1); B(2,1)=quad('asin(x).*x',-1,1); B(3,1)=quad('asin(x).*x.^2',-1,1); B(4,1)=quad('asin(x).*x.^3',-1,1) B

北理工数值分析大作业

数值分析上机作业

第 1 章

1.1计算积分，n=9。（要求计算结果具有6位有效数字）

程序：

n=1:19;

I=zeros(1,19);

I(19)=1/2*((exp(-1)/20)+(1/20));

I(18)=1/2*((exp(-1)/19)+(1/19));

for i=2:10

I(19-i)=1/(20-i)*(1-I(20-i));

end

format long

disp(I(1:19))

结果截图及分析：在MATLAB中运行以上代码，得到结果如下图所示：当计算到数列的第10项时，所得的结果即为n=9时的准确积分值。取6位有效数字可得.

1.2分别将区间[-10.10]分为100,200,400等份，利用mesh或surf

命令画出二元函数

的三维图形。

程序：

>> x = -10:0.1:10;

y = -10:0.1:10;

[X,Y] = meshgrid(x,y);

Z = exp(-abs(X))+cos(X+Y)+1./(X.^2+Y.^2+1);

subplot(2,2,1);

mesh(X,Y,Z);

title('步长0.1')

>> x = -10:0.2:10;

y = -10:0.2:10;

[X,Y] = meshgrid(x,y);

Z = exp(-abs(X))+cos(X+Y)+1./(X.^2+Y.^2+1);

subplot(2,2,1);

mesh(X,Y,Z);

title('步长 0.2')

>>x = -10:0.05:10;

y = -10:0.05:10;

《数值分析》上机实验报告

数值分析上机实验报告

《数值分析》上机实验报告

1.用Newton 法求方程 X 7-X 4+14=0

在（0.1，1.9）中的近似根（初始近似值取为区间端点，迭代6次或误差小于0.00001）。 1.1 理论依据：

设函数在有限区间[a ，b]上二阶导数存在，且满足条件

{}αϕ上的惟一解在区间平方收敛于方程所生的迭代序列

迭代过程由则对任意初始近似值达到的一个中使是其中上不变号

在区间],[0)(3,2,1,0,)

(')

()(],,[x |))(),((|,|,)(||)(|.4;0)(.3],[)(.20

)()(.110......b a x f x k x f x f x x x Newton b a b f a f mir b a c x f a

b c f x f b a x f b f x f k k k k k k ==-

==∈≤-≠>+

令

)9.1()9.1(0)8(4233642)(0)16(71127)(0)9.1(,0)1.0(,1428)(3

333647>⋅''<-=-=''<-=-='<>+-=f f x x x x x f x x x x x f f f x x x f

故以1.9为起点

⎪⎩

⎪⎨

⎧

='-

=+9.1)()(01x x f x f x x k k k k 如此一次一次的迭代，逼近x 的真实根。当前后两个的差<=ε时，就认为求出了近似的根。本程序用Newton 法求代数方程(最高次数不大于10)在（a,b ）区间的根。

1.2 C语言程序原代码：

#include

数值分析上机试题对应参考答案

一、问答题

1、什么是近似值x * 有效数字？

若近似值x*的误差限是某一位的半个单位，该位到x*的第一位非零数字共有n 位，就说x*有n 位有效数字。它可表示为

X=±10m ×(a 1+a 2×10-1+…+a n ×10-(n-1),其中a i (i=1,2,…,n)是0到9中的一个

数字，a 1≠0,m 为整数，且︱x －x *︱≠2

×10m-n+1

2、数值计算应该避免采用不稳定的算法，防止有效数字的损失. 因此，在进行数值运算算法设计过程中主要注意什么？（1）简化计算过程，减少运算次数；（2）避免两个相近的数相减；

（3）避免除数的绝对值远小于被除数的绝对值；（4）防止大数“吃掉”小数的现象；

（5）使用数值稳定的算法，设法控制误差的传播。

3、写出“n 阶阵A 具有n 个不相等的特征值”的等价条件（至少写3 个）

(1)|A|不为零

(2)n 阶矩阵A 的列或行向量组线性无关 (3)矩阵A 为满秩矩阵

(4)n 阶矩阵A 与n 阶可逆矩阵B 等价

4、迭代法的基本思想是什么？

就是用某种极限过程去逐步逼近线性方程组精确解得方法。其基本思想为：先任取一组近似解初值X 0，然后按照某种迭代原则，由X 0计算新的近似解X 1，以此类推，可计算出X 2,X 3，…X K ，。。。,如果｛X ｝收敛，则取为原方程组的解。

5、病态线性方程组的主要判断方法有哪些？

（1）系数矩阵的某两行（列）几乎近似相关（2）系数矩阵的行列式的值很小

（3）用主元消去法解线性方程组时出现小主元

（4）近似解x*已使残差向量r=b-Ax*的范数很小，但该近似解仍不符合问题要求。

数值分析实验

解析解： y(x) x 1 x
附加题
BezΒιβλιοθήκη Baiduer曲线绘制程序设计：
要求：（1）用户输入4个点的坐标，程序根据用户的输入绘制一条Bezier曲线；
（2）增加取样的点数或改变样点的坐标，绘制多条曲线观察数据变化引起的曲线的变化。
（3）当取样点数较多时，采用分段处理的方法绘制曲线（注意连接点处的光滑性）。
要求：
（1）选取不同的初始向量x(0)，给定迭代误差要求
x(k1) x(k ) 105
用Jacobi迭代和Gauss-Seidel迭代法求解，观察得到的序列是否收敛？若收敛，输出结果及迭代次数；
（2）用SOR求解，ω 取1<ω<2中不同值，在
x(k1) x(k ) 105
3 2.1 6
2

x2

5.9
5

2
1 1
5 0
1
2

x3 x4

5

1

要求：用LU分解法和Gauss列主元消元法求解上述方程组，输出Ax=b中A及b，A=LU 中的L及U，detA及向量x；
4、研究解线性方程组Ax=b迭代法收敛速度，给定A∈R20×20为五对角矩阵。
数值分析实验

数值分析上机题

数值剖析上机题

习题 1

17．（上机题）舍入偏差与有效数

1 3 1

设

S N

，其精准值为。

2 2 N

N 1

j 2

（1）编制按从大到小的次序

S N

，计算 S N 的通用程序。

1 3

2 1

N 2

（2）编制按从小到大的次序

S N

，计算 S N 的通用程序。

(N L

N 2 1)2 1

（3）按两种次序分别计算 S

2 ，

S 4 ， S 6

，并指出有效位数。（编制程序时用单精度）

10 10

（4）经过本上机题你理解了什么？

按从大到小的次序计算 S N 的通用程序为：按从小到大的次序计算 S N 的通用程序为：

#include<iostream.h> #include<iostream.h> float sum(float N) float sum(float N) {

{

float j,s,sum=0; float j,s,sum=0; for(j=2;j<=N;j++) for(j=N;j>=2;j--) {

{

s=1/(j*j-1); s=1/(j*j-1); sum+=s;

sum+=s;

}

return sum;

}

从大到小的次序的值

从小到大的次序的值

精准值

有效位数

从大到小

从小到大

5 S 102

S 104

S 106

经过本上机题，看出按两种不一样的次序计算的结果是不同样的，

按从大到小的次序计算

的值与精准值有较大的偏差，而按从小到大的次序计算的值与精准值符合。从大到小的次序

计算获得的结果的有效位数少。