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

合集下载

数值分析上机实验报告

数值分析上机实验报告

实验报告一题目: (绪论) 非线性方程求解及误差估计摘要:非线性方程的解析解通常很难给出,因此线性方程的数值解法就尤为重要。

本实验采用两种常见的求解方法二分法、Newton 法和改进的Newton 法。

可以节省计算机的计算时间,还能减小不必要的误差。

前言:(目的和意义)掌握二分法与Newton 法的基本原理、应用以及熟练掌握用MATLAB 求函数积分 数学原理:(1)函数的调用格式:quadl(filename,a,b,tol,trace)其中filename 是调用函数名,a 和b 分别为定积分的下限和上限。

用来控制积分精度。

(2)秦九韶算法: S n =a nS k =xS k+1+a k (k=n-1,n-2,...,0), P n (x)=S 0 程序设计:例1.1 计算积分de x xx110利用MATLAB ,下面给出主程序>>g=inline('x.^10.*exp(x-1)'); %定义一个语句函数g(x)=exp(x^10*exp(x-1)) I=quadl(g,0,1) I =0.0098例1.9 秦九韶算法a 0=3,a k=2a k-1+3,Pn(x)=a n x^n+a n-1x^(n-1)+...+a1x+a0求I1=P100(0.5),I2=P150(13)>>x=input('x=');n=input('n=');a=3;for i=1:na=2*a+3;ends=z;b=(a-3)/2;for m=1:100s=x*s+b;b=(b-3)/2;enddisp(s);>>x=0.5n=100600.0000>>x=3n=1004.7039e+078结果分析和讨论:结论:对于二分法,只要能够保证在给定的区间内有根,使能够收敛的,当时收敛的速度和给定的区间有关,二且总体上来说速度比较慢。

数值分析第二次大作业——编写ln函数实验报告

数值分析第二次大作业——编写ln函数实验报告

为一个整体,将 ln(a)进行 Taylor 级数分解,进而按式计算可以求得 ln(x)的值。由于我们选 取级数的前 n 项和近似 ln(a),则有:
Rn (a
1)

(1)n
(a 1)n1 n 1
…;
|
Rn
(a
1)
||
(a 1)n1 n 1
|

2
自 92 乔晖 2009011414
精度为 20 位,即 h4 1020 。
int * GetNumber() {return number;}
//获得大数中的整数数组
bool GetSgn() const {return sgn;}
//获得大数的符号
bool IsZero();
//判断大数是否为 0
CBigNumber Reverse() {sgn = !sgn;return *this;} //取相反数
数值分析——编写 ln 函数
实验报告
自 92 乔晖 2009011414
自 92 乔晖 2009011414
数值分析——编写 ln 函数实验报告
目录
一、 实验任务 ................................................... 2
二、 编译环境 ................................................... 2
1u
x 1
35
5
自 92 乔晖 2009011414
数值分析——编写 ln 函数实验报告
其前 n 项和,近似可得 ln(x)的值。
基于以上原理,可以设计迭代算法,累加 2u2i1 得到前 n 项和。每次计算后,判断计 2i 1

数值分析上机实习报告

数值分析上机实习报告

指导教师:姓名:学号:专业:联系电话:上海交通大学目录序言 (3)实验课题(一) 雅可比迭代法和高斯-塞得尔迭代法的收敛性和收敛速度 (4)数值分析 (6)实验课题(二) 松弛因子对SOR法收敛速度的影响 (6)数值分析 (12)总结 (13)附录(程序清单) (14)1.雅可比迭代法和高斯-塞得尔迭代法的收敛性和收敛速度 (14)雅可比迭代法: (14)高斯-塞得尔迭代法: (16)2.松弛因子对SOR法收敛速度的影响 (18)松弛法(SOR) (18)序言随着科学技术的发展,提出了大量复杂的数值计算问题,在实际解决这些计算问题的长期过程中,形成了计算方法这门学科,专门研究各种数学问题的数值解法(近似解法),包括方法的构造和求解过程的误差分析,是一门内容丰富,有自身理论体系的实用性很强的学科。

解决工程问题,往往需要处理很多数学模型,这就要花费大量的人力和时间,但是还有不少数学模型无法用解析法得到解。

使用数值方法并利用计算机,就可以克服这些困难。

事实上,科学计算已经与理论分析、科学实验成为平行的研究和解决科技问题的科学手段,经常被科技工作者所采用。

作为科学计算的核心内容——数值分析(数值计算方法),已逐渐成为广大科技工作者必备的基本知识并越来越被人重视。

由于数值方法是解数值问题的系列计算公式,所以数值方法是否有效,不但与方法本身的好坏有关,而且与数值问题本身的好坏也有关,因此,研究数值方法时,不但需要研究数值方法的好坏,即数值稳定性问题,而且还需要研究数值问题本身的好坏,即数值问题的性态,以及它们的判别问题。

数值计算的绝大部分方法都具有近似性,而其理论又具有严密的科学性,方法的近似值正是建立在理论的严密性基础上,根据计算方法的这一特点。

因此不仅要求掌握和使用算法,还要重视必要的误差分析,以保证计算结果的可靠性。

数值计算还具有应用性强的特点,计算方法的绝大部分方法如求微分方程近似解,求积分近似值,求解超越方程,解线性方程组等都具有较强的实用性,而插值法,最小二乘法,样条函数等也都是工程技术领域中常用的,有实际应用价值的方法。

数值分析实验报告二

数值分析实验报告二

数值实验报告二一、实验名称解线性方程组的列主元素高斯消去法和LU 分解法二、实验目的通过数值实验,从中体会解线性方程组选主元的必要性和LU 分解法的优点,以及方程组系数矩阵和右端向量的微小变化对解向量的影响。

三、实验内容解下列两个线性方程组(1) ⎪⎪⎪⎭⎫ ⎝⎛=⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎭⎫ ⎝⎛--11134.981.4987.023.116.427.199.103.601.3321x x x (2) ⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=⎪⎪⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎪⎪⎭⎫ ⎝⎛----15900001.582012151526099999.23107104321x x x x 四、算法描述1、 列主元素高斯消去法记: ij ij a a =1)( (i, j = 1,2,3n )i i b b =1)( (i = 1,2,3n )消元过程:对于k = 1,2,3n(1) 选行号k i ,使)()(max k i ni k k k i k k a a ≤≤=。

(2) 交换)(k kj a 与)(k j i k a (j = k, k+1,k+2n )以及)()(k i k k k b b 与所含的数值。

(3)对于i = k, k+1,k+2n ,计算)()(k kkk ik ik a a m =)()()1(k kj ik k ij k ij a m a a -=+ (j = k, k+1,k+2n ))()()1(k k ik k i k i b m b b -=+回代过程:)(n nnn n a b x = )()1)()(/(k kk j n k j k kj k k k a x a a x ∑+=-= (k = n-1, n-2, n-3 1 )在此算法中的)(k k i k a 称为第k 个列主元素,它的数值总要被交换到第k 个主对角线元素的位置上。

2、 LU 分解法通过MATLAB 自有的函数,把系数矩阵A 分解成A=LU ,其中:L 是下三角矩阵,U 是上三角矩阵,这时方程组Ax=b 就可以分解成两个容易求解的三角形方程组Ly=b ,Ux=y 。

优质文档精选——数值分析上机实验报告

优质文档精选——数值分析上机实验报告

数值分析上机实验报告《数值分析》上机实验报告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 ab 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)(3225333647>⋅''<-=-=''<-=-='<>+-=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<stdio.h>#include<math.h> main(){double x2,f,f1;double x1=1.9; //取初值为 1.9 do{x2=x1;f=pow(x2,7)-28*pow(x2,4)+14; f1=7*pow(x2,6)-4*28*pow(x2,3); x1=x2-f/f1;}while(fabs(x1-x2)>=0.00001||x1<0.1); //限制循环次数 printf("计算结果:x=%f\n",x1);}1.3 运行结果:1.4 MATLAB 上机程序function y=Newton(f,df,x0,eps,M) d=0;for k=1:Mif feval(df,x0)==0d=2;breakelsex1=x0-feval(f,x0)/feval(df,x0);ende=abs(x1-x0);x0=x1;if e<=eps&&abs(feval(f,x1))<=epsd=1;breakendendif d==1y=x1;elseif d==0y='迭代M次失败';elsey= '奇异'endfunction y=df(x)y=7*x^6-28*4*x^3;Endfunction y=f(x)y=x^7-28*x^4+14;End>> x0=1.9;>> eps=0.00001;>> M=100;>> x=Newton('f','df',x0,eps,M);>> vpa(x,7)1.5 问题讨论:1.使用此方法求方解,用误差来控制循环迭代次数,可以在误差允许的范围内得到比较理想的计算结果。

数值分析上机实验报告

数值分析上机实验报告

数值分析上机实验理学院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 xe 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 εεε+≈22122121+=x x x x x x x )()()(εεε (f(x))()(x)f x εε'≈四.实验内容1. 计算I n=e 1-⎰10nxe x 2dx (n=0,1,...)并估计误差。

解: >> I0 = exp(-1)*quad('(x.^0).*exp(x.^2)',0,1,10^(-10));>> vpa(I0,10) ans =.5380795069>> I1= exp(-1)*quad('(x.^1).*exp(x.^2)',0,1,10^(-10)); >> vpa(I1,10) ans =.3160602794>> I2 = exp(-1)*quad('(x.^2).*exp(x.^2)',0,1,10^(-10)); >> vpa(I2,10) ans =.2309602465>> I3 = exp(-1)*quad('(x.^3).*exp(x.^2)',0,1,10^(-10)); >> vpa(I3,10) ans =.1839397206>> I4 = exp(-1)*quad('(x.^4).*exp(x.^2)',0,1,10^(-10)); >> vpa(I4,10) ans =.1535596302>> I5 = exp(-1)*quad('(x.^5).*exp(x.^2)',0,1,10^(-10)); >> vpa(I5,10) ans =.1321205588>> I6 = exp(-1)*quad('(x.^6).*exp(x.^2)',0,1,10^(-10)); >> vpa(I6,10) ans =.1161009245>> I7 = exp(-1)*quad('(x.^7).*exp(x.^2)',0,1,10^(-10)); >> vpa(I7,10) ans =.1036383235>> I8 = exp(-1)*quad('(x.^8).*exp(x.^2)',0,1,10^(-10)); >> vpa(I8,10) ans =.9364676413e-1>> I9 = exp(-1)*quad('(x.^9).*exp(x.^2)',0,1,10^(-10)); >> vpa(I9,10) ans =.8544670595e-1 2.计算x255的值。

数值分析上机作业2

数值分析上机作业2

数值实验数值实验综述:线性代数方程组的解法是一切科学计算的基础与核心问题。

求解方法大致可分为直接法和迭代法两大类。

直接法——指在没有舍入误差的情况下经过有限次运算可求得方程组的精确解的方法,因此也称为精确法。

当系数矩阵是方的、稠密的、无任何特殊结构的中小规模线性方程组时,Gauss消去法是目前最基本和常用的方法。

如若系数矩阵具有某种特殊形式,则为了尽可能地减少计算量与存储量,需采用其他专门的方法来求解。

Gauss消去等同于矩阵的三角分解,但它存在潜在的不稳定性,故需要选主元素。

对正定对称矩阵,采用平方根方法无需选主元。

方程组的性态与方程组的条件数有关,对于病态的方程组必须采用特殊的方法进行求解。

实验一一、实验名称:矩阵的LU分解.二、实验目的:用不选主元的LU分解和列主元LU分解求解线性方程组Ax=b, 并比较这两种方法.三、实验内容与要求(1)用所熟悉的计算机语言将不选主元和列主元LU分解编成通用的子程序,然后用编写的程序求解下面的84阶方程组将计算结果与方程组的精确解进行比较,并就此谈谈你对Gauss消去法的看法。

(2)写出追赶法求解三对角方程组的过程,并编写程序求该实验中的方程组(1)①不选主元高斯消去法求解方程组function [L,U]=gauss1(A,b)n=length(A);for k=1:n-1A(k+1:n,k)=A(k+1:n,k)/A(k,k);A(k+1:n,k+1:n)=A(k+1:n,k+1:n)-A(k+1:n,k)*A(k,k+1:n);endL=tril(A(:,1:n),-1)+eye(n);U=triu(A);主程序为:Clear;clc;a=ones([84,1])*6;b=ones([83,1]);c=ones([83,1])*8;A=diag(a)+diag(b,1)+diag(c,-1); d=ones([82,1])*15;b=[7;d;14];[L U]=gauss1(A,b);n=length(A);y=ones(n,1);y=L\b;x=ones(n,1);x=U\y解为:x=11.000000000000001.000000000000001.000000000000000.9999999999999981.000000000000000.9999999999999931.000000000000010.9999999999999721.000000000000060.9999999999998861.000000000000230.9999999999995451.000000000000910.9999999999981811.000000000003640.9999999999927251.000000000014550.9999999999708981.000000000058200.9999999998835921.000000000232820.9999999995343671.000000000931270.9999999981374691.000000003725060.9999999925498741.000000014900250.9999999701994971.000000059601010.9999998807979861.000000238404030.9999995231919461.000000953616110.9999980927677831.000003814464440.9999923710711301.000015257857740.9999694842845201.000061031430960.9998779371380811.000244125723840.9995117485523231.000976502895350.9980469942092911.003906011581410.9921879768371871.015624046325570.9687519073490881.062496185300920.8750076294018071.249984741181830.5000305176945351.99993896437811 -0.999877927824961 4.99975585192486 -6.99951168894947 16.9990233182979 -30.9980463981918 64.9960918427676 -126.992179871071 256.984344484284 -510.968627937136 1024.93701174855 -2046.87304699420 4096.74218797682 -8190.46875190732 16383.8750076293 -32764.5000305175 65531.0001220701 -131055.000488281 262097.001953124 -524127.007812498 1048001.03125000 -2094975.124999994185857.49999999-8355328.9999999716645128.9999999-33028126.999999965007744.9999998-125821439.000000234866688.999999-402628606.999999536838144.999998可见,这是一个病态方程,从56个跟开始发散。

数值分析上机实践报告

数值分析上机实践报告

数值分析上机实践报告一、实验目的本次实验主要目的是通过上机操作,加深对数值分析算法的理解,并熟悉使用Matlab进行数值计算的基本方法。

在具体实验中,我们将实现三种常见的数值分析算法:二分法、牛顿法和追赶法,分别应用于解决非线性方程、方程组和线性方程组的求解问题。

二、实验原理与方法1.二分法二分法是一种常见的求解非线性方程的数值方法。

根据函数在给定区间端点处的函数值的符号,不断缩小区间的长度,直到满足精度要求。

2.牛顿法牛顿法是求解方程的一种迭代方法,通过构造方程的泰勒展开式进行近似求解。

根据泰勒展式可以得到迭代公式,利用迭代公式不断逼近方程的解。

3.追赶法追赶法是用于求解三对角线性方程组的一种直接求解方法。

通过构造追赶矩阵,采用较为简便的向前追赶和向后追赶的方法进行计算。

本次实验中,我们选择了一组非线性方程、方程组和线性方程组进行求解。

具体的实验步骤如下:1.调用二分法函数,通过输入给定区间的上下界、截止误差和最大迭代次数,得到非线性方程的数值解。

2.调用牛顿法函数,通过输入初始迭代点、截止误差和最大迭代次数,得到方程组的数值解。

3.调用追赶法函数,通过输入追赶矩阵的三个向量与结果向量,得到线性方程组的数值解。

三、实验结果与分析在进行实验过程中,我们分别给定了不同的参数,通过调用相应的函数得到了实验结果。

下面是实验结果的汇总及分析。

1.非线性方程的数值解我们通过使用二分法对非线性方程进行求解,给定了区间的上下界、截止误差和最大迭代次数。

实验结果显示,根据给定的输入,我们得到了方程的数值解。

通过与解析解进行比较,可以发现二分法得到的数值解与解析解的误差在可接受范围内,说明二分法是有效的。

2.方程组的数值解我们通过使用牛顿法对方程组进行求解,给定了初始迭代点、截止误差和最大迭代次数。

实验结果显示,根据给定的输入,我们得到了方程组的数值解。

与解析解进行比较,同样可以发现牛顿法得到的数值解与解析解的误差在可接受范围内,说明牛顿法是有效的。

哈工大数值分析报告上机实验报告

哈工大数值分析报告上机实验报告

实验报告一题目: Gauss 列主元消去法摘要:求解线性方程组地方法很多,主要分为直接法和间接法.本实验运用直接法地Guass 消去法,并采用选主元地方法对方程组进行求解.前言:(目地和意义)1. 学习Gauss 消去法地原理.2. 了解列主元地意义.3. 确定什么时候系数阵要选主元数学原理:由于一般线性方程在使用Gauss 消去法求解时,从求解地过程中可以看到,若)1(-k kk a =0,则必须进行行交换,才能使消去过程进行下去.有地时候即使≠-)1(k kk a 0,但是其绝对值非常小,由于机器舍入误差地影响,消去过程也会出现不稳定得现象,导致结果不正确.因此有必要进行列主元技术,以最大可能地消除这种现象.这一技术要寻找行r ,使得)1()1(max ||->-=k ik ki k rk a a 并将第r 行和第k 行地元素进行交换,以使得当前地)1(-k kk a 地数值比0要大地多.这种列主元地消去法地主要步骤如下:1. 消元过程对k =1,2,…,n -1,进行如下步骤.1) 选主元,记ik ki rk a a >=max || 若||rk a 很小,这说明方程地系数矩阵严重病态,给出警告,提示结果可能不对.2) 交换增广阵A 地r ,k 两行地元素.kj rj a a ↔ (j=k,…,n +1)3) 计算消元kk kj ik ij ij a a a a a /-= (i=k+1,…,n ; j =k +1,……,n +1)2. 回代过程对k = n , n -1,…,1,进行如下计算)/(11,∑-=+-=nk j kk j kj n k k a x a a x至此,完成了整个方程组地求解.程序设计:本实验采用Matlab地M文件编写.Gauss消去法源程序:cleara=input('输入系数阵:>>\n')b=input('输入列阵b:>>\n')n=length(b);A=[a b]x=zeros(n,1);%%%函数主体for k=1:n-1;%%%是否进行主元选取if abs(A(k,k))<yipusilong;%事先给定地认为有必要选主元地小数yzhuyuan=1;else yzhuyuan=0;endif yzhuyuan;%%%%选主元t=A(k,k);for r=k+1:n;if abs(A(r,k))>abs(t)p=r;else p=k;endend%%%交换元素if p~=k;for q=k:n+1;s=A(k,q);A(k,q)=A(p,q);A(p,q)=s;endendend%%%判断系数矩阵是否奇异或病态非常严重if abs(A(k,k))< yipusilongdisp(‘矩阵奇异,解可能不正确’)end%%%%计算消元,得三角阵for r=k+1:n;m=A(r,k)/A(k,k);for q=k:n+1;A(r,q)=A(r,q)-A(k,q)*m;endendend%%%%求解xx(n)=A(n,n+1)/A(n,n);for k=n-1:-1:1;s=0;for r=k+1:n;s=s+A(k,r)*x(r);endt=(A(k,n+1)-s)x(k)=(A(k,n+1)-s)/A(k,k)end结果分析和讨论:例:求解方程⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡321643.5072.12-623.4712.31-32108-z y x .求解地结果为:x =[]367257386.0,05088607.0-49105822.0-, 例:求解方程⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡73109104-10172-42-4z y x 求得地结果为:x =[]857142857.1,89285714.0-196428571.0, 结论:采用Gauss 消去法时,如果在消元时对角线上地元素始终较大(假如大于10-5),那么本方法不需要进行列主元计算,计算结果一般就可以达到要求,否则必须进行列主元这一步,以减少机器误差带来地影响,使方法得出地结果正确.实验报告二题目: Rung 现象产生和克服摘要:由于高次多项式插值不收敛,会产生Runge 现象,本实验在给出具体地实例后,采用分段线性插值和三次样条插值地方法有效地克服了这一现象,而且还取地很好地插值效果.前言:(目地和意义)1. 深刻认识多项式插值地缺点.2. 明确插值地不收敛性怎样克服.3. 明确精度与节点和插值方法地关系.数学原理:在给定n+1个节点和相应地函数值以后构造n 次地Lagrange 插值多项式,实验结果表明(见后面地图)这种多项式并不是随着次数地升高对函数地逼近越来越好,这种现象就是Rung 现象.解决Rung 现象地方法通常有分段线性插值、三次样条插值等方法.分段线性插值:设在区间[a, b ]上,给定n+1个插值节点a=x 0<x 1<…<x n =b和相应地函数值y 0,y 1,…,y n ,,求作一个插值函数)(x φ,具有如下性质:1) j j y x =)(φ,j=0,1,…,n .2) )(x φ在每个区间[x i , x j ]上是线性连续函数.则插值函数)(x φ称为区间[a, b ]上对应n 个数据点地分段线性插值函数.三次样条插值:给定区间[a, b ]一个分划⊿:a=x 0<x 1<…<x N =b若函数S(x)满足下列条件:1) S(x)在每个区间[x i , x j ]上是不高于3次地多项式.2) S(x)及其2阶导数在[a, b ]上连续.则称S(x)使关于分划⊿地三次样条函数. 程序设计流程:本实验采用Matlab 地M 文件编写.其中待插值地方程写成function 地方式,如下function y=f(x);y=1/(1+25*x*x );写成如上形式即可,下面给出主程序Lagrange 插值源程序:n=input('将区间分为地等份数输入:\n');s=[-1+2/n*[0:n]];%%%给定地定点,Rf为给定地函数x=-1:0.01:1;f=0;for q=1:n+1;l=1;%求插值基函数for k=1:n+1;if k~=q;l=l.*(x-s(k))./(s(q)-s(k));elsel=l;endendf=f+Rf(s(q))*l;%求插值函数endplot(x,f,'r')%作出插值函数曲线grid onhold on分段线性插值源程序clearn=input('将区间分为地等份数输入:\n');s=[-1+2/n*[0:n]];%%%给定地定点,Rf为给定地函数m=0;hh=0.001;for x=-1:hh:1;ff=0;for k=1:n+1;%%%求插值基函数switch kcase 1if x<=s(2);l=(x-s(2))./(s(1)-s(2));elsel=0;endcase n+1if x>s(n);l=(x-s(n))./(s(n+1)-s(n));elsel=0;endotherwiseif x>=s(k-1)&x<=s(k);l=(x-s(k-1))./(s(k)-s(k-1));else if x>=s(k)&x<=s(k+1);l=(x-s(k+1))./(s(k)-s(k+1));elsel=0;endendendff=ff+Rf(s(k))*l;%%求插值函数值endm=m+1;f(m)=ff;end%%%作出曲线x=-1:hh:1;plot(x,f,'r');grid onhold on三次样条插值源程序:(采用第一边界条件)clearn=input('将区间分为地等份数输入:\n');%%%插值区间a=-1;b=1;hh=0.001;%画图地步长s=[a+(b-a)/n*[0:n]];%%%给定地定点,Rf为给定地函数%%%%第一边界条件Rf"(-1),Rf"(1)v=5000*1/(1+25*a*a)^3-50/(1+25*a*a)^4;for k=1:n;%取出节点间距h(k)=s(k+1)-s(k);endfor k=1:n-1;%求出系数向量lamuda,miula(k)=h(k+1)/(h(k+1)+h(k));miu(k)=1-la(k);end%%%%赋值系数矩阵Afor k=1:n-1;for p=1:n-1;switch pcase kA(k,p)=2;case k-1A(k,p)=miu(p+1);case k+1A(k,p)=la(p-1);otherwiseA(k,p)=0;endendend%%%%求出d阵for k=1:n-1;switch kcase 1d(k)=6*f2c([s(k) s(k+1) s(k+2)])-miu(k)*v;case n-1d(k)=6*f2c([s(k) s(k+1) s(k+2)])-la(k)*v;otherwised(k)=6*f2c([s(k) s(k+1) s(k+2)]);endend%%%%求解M阵M=A\d';M=[v;M;v];%%%%m=0;f=0;for x=a:hh:b;if x==a;p=1;elsep=ceil((x-s(1))/((b-a)/n));endff1=0;ff2=0;ff3=0;ff4=0;m=m+1;ff1=1/h(p)*(s(p+1)-x)^3*M(p)/6;ff2=1/h(p)*(x-s(p))^3*M(p+1)/6;ff3=((Rf(s(p+1))-Rf(s(p)))/h(p)-h(p)*(M(p+1)-M(p))/6)*(x-s(p));ff4=Rf(s(p))-M(p)*h(p)*h(p)/6;f(m)=ff1+ff2+ff3+ff4 ;end%%%作出插值图形x=a:hh:b;plot(x,f,'k')hold ongrid on结果分析和讨论: 本实验采用函数22511)(xx f +=进行数值插值,插值区间为[-1,1],给定节点为 x j =-1+jh ,h=0.1,j =0,…,n .下面分别给出Lagrang e 插值,三次样条插值,线性插值地函数曲线和数据表.图中只标出Lagrang e 插值地十次多项式地曲线,其它曲线没有标出,从数据表中可以看出具体地误差.表中,L10(x)为Lagrang e插值地10次多项式,S10(x),S40(x)分别代表n=10,40地三次样条插值函数,X10(x),X40(x)分别代表n=10,40地线性分段插值函数.x f(x)L10(x)S10(x) S40(x) X10(x) X40(x) -1.00000000000000 0.03846153846154 0.03846153846154 0.03846153846154 0.03846153846154 0.03846153846154 0.03846153846154 -0.95000000000000 0.04244031830239 1.92363114971920 0.04240833151040 0.04244031830239 0.04355203619910 0.04244031830239 -0.90000000000000 0.04705882352941 1.57872099034926 0.04709697585458 0.04705882352941 0.04864253393665 0.04705882352941 -0.85000000000000 0.05245901639344 0.71945912837982 0.05255839923979 0.05245901639344 0.05373303167421 0.05245901639344 -0.80000000000000 0.05882352941176 0.05882352941176 0.05882352941176 0.05882352941176 0.05882352941176 0.05882352941176 -0.75000000000000 0.06639004149378 -0.23146174989674 0.06603986172744 0.06639004149378 0.06911764705882 0.06639004149378 -0.70000000000000 0.07547169811321 -0.22619628906250 0.07482116198866 0.07547169811321 0.07941176470588 0.07547169811321 -0.65000000000000 0.08648648648649 -0.07260420322418 0.08589776360849 0.08648648648649 0.08970588235294 0.08648648648649 -0.60000000000000 0.10000000000000 0.10000000000000 0.10000000000000 0.10000000000000 0.10000000000000 0.10000000000000 -0.55000000000000 0.11678832116788 0.21559187891257 0.11783833017713 0.11678832116788 0.12500000000000 0.11678832116788 -0.50000000000000 0.13793103448276 0.25375545726103 0.14004371555730 0.13793103448276 0.15000000000000 0.13793103448276 -0.45000000000000 0.16494845360825 0.23496854305267 0.16722724315883 0.16494845360825 0.17500000000000 0.16494845360825 -0.40000000000000 0.20000000000000 0.20000000000000 0.20000000000000 0.20000000000000 0.20000000000000 0.20000000000000 -0.35000000000000 0.24615384615385 0.19058046675376 0.24054799403464 0.24615384615385 0.27500000000000 0.24615384615385 -0.30000000000000 0.30769230769231 0.23534659131080 0.29735691695860 0.30769230769231 0.35000000000000 0.30769230769231 -0.25000000000000 0.39024390243902 0.34264123439789 0.38048738140327 0.39024390243902 0.42500000000000 0.39024390243902 -0.20000000000000 0.50000000000000 0.50000000000000 0.50000000000000 0.50000000000000 0.50000000000000 0.50000000000000 -0.15000000000000 0.64000000000000 0.67898957729340 0.65746969368431 0.64000000000000 0.62500000000000 0.64000000000000 -0.10000000000000 0.80000000000000 0.84340742982890 0.82052861660828 0.80000000000000 0.75000000000000 0.80000000000000 -0.05000000000000 0.94117647058824 0.95862704866073 0.94832323122810 0.94117647058824 0.87500000000000 0.941176470588240 1.00000000000000 1.00000000000000 1.00000000000000 1.00000000000000 1.000000000000001.000000000000000.05000000000000 0.94117647058824 0.95862704866073 0.94832323122810 0.94117647058824 0.87500000000000 0.941176470588240.10000000000000 0.80000000000000 0.84340742982890 0.82052861660828 0.80000000000000 0.75000000000000 0.800000000000000.15000000000000 0.64000000000000 0.67898957729340 0.65746969368431 0.64000000000000 0.62500000000000 0.640000000000000.20000000000000 0.50000000000000 0.50000000000000 0.50000000000000 0.50000000000000 0.50000000000000 0.500000000000000.25000000000000 0.39024390243902 0.34264123439789 0.38048738140327 0.39024390243902 0.42500000000000 0.390243902439020.30000000000000 0.30769230769231 0.23534659131080 0.29735691695860 0.30769230769231 0.350000000000000.307692307692310.35000000000000 0.24615384615385 0.19058046675376 0.24054799403464 0.24615384615385 0.27500000000000 0.246153846153850.40000000000000 0.20000000000000 0.20000000000000 0.20000000000000 0.20000000000000 0.20000000000000 0.200000000000000.45000000000000 0.16494845360825 0.23496854305267 0.16722724315883 0.16494845360825 0.17500000000000 0.164948453608250.50000000000000 0.13793103448276 0.25375545726103 0.14004371555730 0.13793103448276 0.15000000000000 0.137931034482760.55000000000000 0.11678832116788 0.21559187891257 0.11783833017713 0.11678832116788 0.12500000000000 0.116788321167880.60000000000000 0.10000000000000 0.10000000000000 0.10000000000000 0.10000000000000 0.10000000000000 0.100000000000000.65000000000000 0.08648648648649 -0.07260420322418 0.08589776360849 0.08648648648649 0.08970588235294 0.086486486486490.70000000000000 0.07547169811321 -0.22619628906250 0.07482116198866 0.07547169811321 0.07941176470588 0.075471698113210.75000000000000 0.06639004149378 -0.23146174989674 0.06603986172744 0.06639004149378 0.06911764705882 0.066390041493780.80000000000000 0.05882352941176 0.05882352941176 0.05882352941176 0.05882352941176 0.05882352941176 0.058823529411760.85000000000000 0.05245901639344 0.71945912837982 0.05255839923979 0.05245901639344 0.05373303167421 0.052459016393440.90000000000000 0.04705882352941 1.57872099034926 0.04709697585458 0.04705882352941 0.04864253393665 0.047058823529410.95000000000000 0.04244031830239 1.92363114971920 0.04240833151040 0.04244031830239 0.04355203619910 0.042440318302391.00000000000000 0.03846153846154 0.03846153846154 0.03846153846154 0.03846153846154 0.03846153846154 0.03846153846154从以上结果可以看到,用三次样条插值和线性分段插值,不会出现多项式插值是出现地Runge现象,插值效果明显提高.进一步说,为了提高插值精度,用三次样条插值和线性分段插值是可以增加插值节点地办法来满足要求,而用多项式插值函数时,节点数地增加必然会使多项式地次数增加,这样会引起数值不稳定,所以说这两种插值要比多项式插值好地多.而且在给定节点数地条件下,三次样条插值地精度要优于线性分段插值,曲线地光滑性也要好一些.实验报告三题目: 多项式最小二乘法摘要:对于具体实验时,通常不是先给出函数地解析式,再进行实验,而是通过实验地观察和测量给出离散地一些点,再来求出具体地函数解析式.又因为测量误差地存在,实际真实地解析式曲线并不一定通过测量给出地所有点.最小二乘法是求解这一问题地很好地方法,本实验运用这一方法实现对给定数据地拟合. 前言:(目地和意义)1. 学习使用最小二成法地原理2. 了解法方程地特性 数学原理:对于给定地测量数据(x i ,f i )(i=1,2,…,n ),设函数分布为∑==mj j j x a x y 0)()(ϕ特别地,取)(x j ϕ为多项式j j x x =)(ϕ (j=0, 1,…,m )则根据最小二乘法原理,可以构造泛函∑∑==-=n i mj i j j i m x a f a a a H 110))((),,,(ϕ令0=∂∂ka H(k=0, 1,…,m ) 则可以得到法方程⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡),(),(),(),(),(),(),(),(),(),(),(),(1010101111000100m m m m m m m m f f f a a a ϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕϕ求该解方程组,则可以得到解m a a a ,,,10 ,因此可得到数据地最小二乘解∑=≈mj j j x a x f 0)()(ϕ程序设计:本实验采用Matlab 地M 文件编写.其中多项式函数j j x =ϕ写成function 地方式,如下function y=fai(x,j)y=1;for i=1:jy=x.*y;end写成如上形式即可,下面给出主程序.多项式最小二乘法源程序clear%%%给定测量数据点(s,f)s=[3 4 5 6 7 8 9];f=[2.01 2.98 3.50 5.02 5.47 6.02 7.05];%%%计算给定地数据点地数目n=length(f);%%%给定需要拟合地数据地最高次多项式地次数m=10;%%%程序主体for k=0:m;g=zeros(1,m+1);for j=0:m;t=0;for i=1:n;%计算内积(fai(si),fai(si))t=t+fai(s(i),j)*fai(s(i),k);endg(j+1)=t;endA(k+1,:)=g;%法方程地系数矩阵t=0;for i=1:n;%计算内积(f(si),fai(si))t=t+f(i)*fai(s(i),k);endb(k+1,1)=t;enda=A\b%求出多项式系数x=[s(1):0.01:s(n)]';y=0;for i=0:m;y=y+a(i+1)*fai(x,i);endplot(x,y)%作出拟合成地多项式地曲线grid onhold onplot(s,f,'rx') %在上图中标记给定地点结果分析和讨论:例用最小二乘法处理下面地实验数据.并作出)f地近似分布图.(x分别采用一次,二次和五次多项式来拟合数据得到相应地拟合多项式为:y1=-0.38643+0.82750x;y2=-1.03024+1.06893x-0.02012x2;y5=-50.75309+51.53527x-19.65947x2+3.66585x3-0.32886x4+0.01137x5;分别作出它们地曲线图,图中点划线为y1曲线,实线为y2曲线,虚线为y5曲线.’x’为给定地数据点.从图中可以看出并不是多项式次数越高越好,次数高了,曲线越能给定点处和实际吻合,但别地地方就很差了.因此,本例选用一次和两次地多项式拟合应该就可以了.实验报告四题目: Romberg 积分法摘要:对于实际地工程积分问题,很难应用Newton-Leibnitz 公式去求解.因此应用数值方法进行求解积分问题已经有着很广泛地应用,本文基于Romberg 积分法来解决一类积分问题.前言:(目地和意义)1. 理解和掌握Romberg 积分法地原理;2. 学会使用Romberg 积分法;3. 明确Romberg 积分法地收敛速度及应用时容易出现地问题. 数学原理:考虑积分⎰=ba dx x f f I )()(,欲求其近似值,通常有复化地梯形公式、Simpsion公式和Cotes 公式.但是给定一个精度,这些公式达到要求地速度很缓慢.如何提高收敛速度,自然是人们极为关心地课题.为此,记T 1,k 为将区间[a,b ]进行2k 等分地复化地梯形公式计算结果,记T 2,k 为将区间[a,b ]进行2k 等分地复化地Simpsion 公式计算结果,记T 3,k 为将区间[a,b ]进行2k 等分地复化地Cotes 公式计算结果.根据Richardson 外推加速方法,可以得到收敛速度较快地Romberg 积分法.其具体地计算公式为: 1. 准备初值,计算)]()([21,1b f a f ba T +-=2. 按梯形公式地递推关系,计算∑-=-+-+-+-+=1201,11,11))5.0(2(221k i k k k k i ab a f a b T T 3. 按Romberg 积分公式计算加速值1441,11,11,--=----+---m mk m m k m m m k m T T T m=2,…,k4. 精度控制.对给定地精度R ,若R T T m m <--1,11,则终止计算,并取1,m T 为所求结果;否则返回2重复计算,直至满足要求地精度为止. 程序设计:本实验采用Matlab 地M 文件编写.其中待积分地函数写成function 地方式,例如如下function yy=f(x,y); yy=x.^3;写成如上形式即可,下面给出主程序Romberg 积分法源程序%%% Romberg 积分法 clear%%%积分区间 b=3; a=1;%%%精度要求 R=1e-5;%%%应用梯形公式准备初值 T(1,1)=(b-a)*(f(b)+f(a))/2; T(1,2)=T(1,1)/2+(b-a)/2*f((b+a)/2); T(2,1)=(4*T(1,2)-T(1,1))/(4-1); j=2; m=2;%%%主程序体%%%while(abs(T(m,1)-T(m-1,1))>R);%%%精度控制 j=j+1; s=0;for p=1:2^(j-2);s=s+f(a+(2*p-1)*h/(2^(j-1))); endT(1,j)=T(1,j-1)/2+h*s/(2^(j-1)); %%%梯形公式应用 for m=2:j; k=(j-m+1);T(m,k)=((4^(m-1))*T(m-1,k+1)-T(m-1,k))/(4^(m-1)-1); end end%%%给出 Romberg 积分法地函数表 I=T(m,1)结果分析和讨论: 1. 求积分dx x10063.精确解I= 24999676.运行程序得Romberg 积分法地函数表为1.0e+007 *4.70101520000000 3.05022950000000 2.63753307500000 2.49996760000000 2.49996760000000 0 2.49996760000000 0 0由函数表知Romberg 积分给出地结果为2.4999676*10^7,与精确没有误差,精度很高.2. 求积分dx xx⎰10sin . 直接按前面方法进行积分,会发现系统报错,出现了0为除数地现象.出现这种情况地原因就是当x=0时,被积函数分母出现0,如果用一个适当地小数ε(最好不要小于程序给定地最小误差值,但不能小于机器地最大精度)来代替可以避免这个问题.本实验取R =ε,可得函数表为:0.92073548319659 0.93979327500190 0.94451351171417 0.94569085359489 0.94598501993743 0.94614587227034 0.94608692395160 0.94608330088846 0.94608307538495 0 0.94608299406368 0.94608305935092 0.94608306035138 0 0 0.94608306038722 0.94608306036726 0 0 0 0.94608306036718 0 0 0 0故该函数地积分为0.94608306036718,取8位有效数字.3. 求积分dx x ⎰12sin本题地解析解很难给出,但运用Romberg 积分可以很容易给出近似解,函数表为:0.42073549240395 0.33406972582924 0.31597536075922 0.31168023948094 0.31062036680949 0.31035626065456 0.30518113697100 0.30994390573588 0.31024853238818 0.31026707591900 0.31026822526959 0 0.31026142365354 0.31026884083167 0.31026831215439 0.31026830189296 0 0 0.31026895856465 0.31026830376269 0.31026830173008 0 0 0 0.31026830119484 0.31026830172211 0 0 0 0 0.31026830172262 0 0 0 0 0故该函数地积分为0.31026830172262,取8位有效数字.结论:Romberg 积分通常要求被积函数在积分区间上没有奇点.如有奇点,且奇点为第一间断点,那么采用例3地方法,还是能够求出来地,否则,必须采用其它地积分方法.当然,Romberg 积分地收敛速度还是比较快地.。

数值分析上机实习报告

数值分析上机实习报告

数值分析上机实习报告随着现代科学技术的迅猛发展,计算机科学的应用日益广泛,数值分析作为计算机科学中重要的分支之一,其在工程、物理、生物学等领域的应用也越来越受到重视。

本学期,我们在数值分析课程的学习中,进行了多次上机实习,通过实习,我们对数值分析的基本方法和算法有了更深入的理解和掌握。

在实习过程中,我们使用了MATLAB软件作为主要的工具,MATLAB是一种功能强大的数学软件,它提供了丰富的数值计算函数和图形显示功能,使我们能够更加方便地进行数值计算和分析。

第一次实习是线性插值和函数逼近。

我们学习了利用已知数据点构造插值函数的方法,并通过MATLAB软件实现了线性插值和拉格朗日插值。

通过实习,我们了解了插值的基本原理,掌握了插值的计算方法,并能够利用MATLAB软件进行插值计算。

第二次实习是解线性方程组。

我们学习了高斯消元法、列主元高斯消元法和克莱姆法则等解线性方程组的方法,并通过MATLAB软件实现了这些算法。

在实习过程中,我们通过实际例子了解了这些算法的应用,掌握了它们的计算步骤,并能够利用MATLAB软件准确地求解线性方程组。

第三次实习是求解非线性方程和方程组。

我们学习了二分法、牛顿法、弦截法和迭代法等求解非线性方程的方法,以及雅可比法和高斯-赛德尔法等求解非线性方程组的方法。

通过实习,我们了解了非线性方程和方程组的求解方法,掌握了它们的计算步骤,并能够利用MATLAB软件求解实际问题。

通过这次上机实习,我们不仅深入学习了数值分析的基本方法和算法,而且锻炼了利用MATLAB软件进行数值计算和分析的能力。

同时,我们也认识到了数值分析在实际问题中的应用价值,增强了解决实际问题的能力。

总之,这次上机实习使我们受益匪浅,对我们学习数值分析课程起到了很好的辅助作用。

数值分析实验报告

数值分析实验报告

《数值分析》实验报告班级:姓名:学号:指导老师:实验基本要求一、上机前的准备工作1、复习和掌握与本次实验有关的教学内容。

2、根据本次实验要求,在纸上编写算法及上机的程序,并经过人工模拟运行检验,减少不必要的错误,提高上机效率。

切忌不编程序、不作人工检查就进行程序输入,这只能使上机调试的难度增加,甚至可能带来学习自信心的下降,影响后续课程的学习。

二、上机实验步骤1、启动开发环境;2、建立源程序文件,输入源程序;3、编译产生目标程序,连接生成可执行程序,运行程序,输出结果;4、对数值计算结果进行误差分析,讨论数值算法的收敛性与稳定性;5、整理实验报告。

三、实验报告实验报告是记录实验工作全过程的技术文档,实验报告的撰写是科学技术工作的一个组成部分。

《数值分析》实验报告包括下列要求:1、实验原理;2、实验内容和要求;3、数值算法描述,包括数据输入、数据处理和数据输出;4、算法的实现(1)给出具体的计算实例,(2)经调试正确的源程序清单,(3)对具体的数值例子给出数值结果;5、计算结果的误差分析,算法的收敛性与稳定性的讨论;6、实验心得。

实验一、误差分析一、实验目的1、通过上机编程,复习巩固以前所学程序设计语言及上机操作指令;2、通过上机计算,了解误差、绝对误差、误差界、相对误差界的有关概念;3、 通过上机计算,了解舍入误差所引起的数值不稳定性。

二、实验原理误差问题是数值分析的基础,又是数值分析中一个困难的课题。

在实际计算中,如果选用了不同的算法,由于舍入误差的影响,将会得到截然不同的结果。

因此,选取算法时注重分析舍入误差的影响,在实际计算中是十分重要的。

同时,由于在数值求解过程中用有限的过程代替无限的过程会产生截断误差,因此算法的好坏会影响到数值结果的精度。

三、实验任务对20,,2,1,0 =n ,计算定积分⎰+=105dx x x y nn . 算法1:利用递推公式151--=n n y ny , 20,,2,1 =n , 取 ⎰≈-=+=100182322.05ln 6ln 51dx x y . 算法2:利用递推公式n n y n y 51511-=- 1,,19,20 =n . 注意到⎰⎰⎰=≤+≤=1010202010201051515611261dx x dx x x dx x , 取 008730.0)12611051(20120≈+≈y . 思考:从计算结果看,哪个算法是不稳定的,哪个算法是稳定的。

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

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

数值分析上机实验报告《数值分析》上机实验报告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 ab 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)(3225333647>⋅''<-=-=''<-=-='<>+-=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<stdio.h>#include<math.h>main(){double x2,f,f1;double x1=1.9; //取初值为1.9do{x2=x1;f=pow(x2,7)-28*pow(x2,4)+14;f1=7*pow(x2,6)-4*28*pow(x2,3);x1=x2-f/f1;}while(fabs(x1-x2)>=0.00001||x1<0.1); //限制循环次数printf("计算结果:x=%f\n",x1);}1.3 运行结果:1.4 MATLAB上机程序function y=Newton(f,df,x0,eps,M)d=0;for k=1:Mif feval(df,x0)==0d=2;breakelsex1=x0-feval(f,x0)/feval(df,x0);ende=abs(x1-x0);x0=x1;if e<=eps&&abs(feval(f,x1))<=epsd=1;breakendendif d==1y=x1;elseif d==0y='迭代M次失败';elsey= '奇异'endfunction y=df(x)y=7*x^6-28*4*x^3;Endfunction y=f(x)y=x^7-28*x^4+14;End>> x0=1.9;>> eps=0.00001;>> M=100;>> x=Newton('f','df',x0,eps,M);>> vpa(x,7)1.5 问题讨论:1.使用此方法求方解,用误差来控制循环迭代次数,可以在误差允许的范围内得到比较理想的计算结果。

《数值分析》实验报告2

《数值分析》实验报告2

《数值分析》实验报告一、问题的提出求解线性方程组的迭代法,即是用某种极限过程去逐步逼近线性方程组的精确解的过程,迭代法是解大型稀疏矩阵方程组的重要方法。

二、实验名称运用MATLAB编程实现雅可比(Jacobi)迭代和高斯-赛德尔(Gauss-Seidel)迭代。

三、实验目的1、熟悉了解雅可比(Jacobi)迭代和高斯-赛德尔(Gauss-Seidel)迭代的算法。

2、学习MATLAB软件的功能。

四、基本原理五、实验环境操作环境:Windows10实验平台:Matlab7.1软件六、试验设计1、jacobi迭代法(1)算例:课本p54页例1(2)程序清单Jacobi迭代法的MATLAB函数文件Jacobi.m如下:function [y,n]=jacobi(A,b,x0,eps)if nargin==3eps=1.0e-6;elseif nargin<3errorreturnendD=diag(diag(A)); %求A的对角矩阵L=-tril(A,-1); %求A的下三角阵U=-triu(A,1); %求A的上三角阵B=D\(L+U);f=D\b;y=B*x0+f;n=1; %迭代次数while norm(y-x0)>=epsx0=y;y=B*x0+f;n=n+1;end(3)实验结果及分析:>> A=[8,-3,2;4,11,-1;6,3,12];>> b=[20,33,36]';>> [x,n]=jacobi(A,b,[0,0,0]',1.0e-6)x =3.00002.00001.0000n =162、Gauss-seidel迭代法(1)算例:课本p54页例1(2)程序清单:Gauss-Serdel迭代法的MATLAB函数文件gauseidel.m如下:function [y,n]=gauseidel(A,b,x0,eps)if nargin==3eps=1.0e-6;elseif nargin<3errorreturnendD=diag(diag(A)); %求A的对角矩阵L=-tril(A,-1); %求A的下三角阵U=-triu(A,1); %求A的上三角阵G=(D-L)\U;f=(D-L)\b;y=G*x0+f;n=1; %迭代次数while norm(y-x0)>=epsx0=y;y=G*x0+f;n=n+1;end(3)实验结果及分析:>> A=[8,-3,2;4,11,-1;6,3,12];>> b=[20,33,36]';>> [x,n]=gauseidel(A,b,[0,0,0]',1.0e-6)x =3.00002.00001.0000n =9七、结果说明:高斯-赛德尔迭代法比雅可比迭代法收敛得快一些(达到相同精度所需迭代次数较少)。

数值分析上机作业(2)

数值分析上机作业(2)

一、数值求解如下正方形域上的Poisson 方程边值问题 2222(,)1,0,1(0,)(1,)(1),01(,0)(,1)0,01u u f x y x y x y u y u y y y y u x u x x ⎧⎛⎫∂∂-+==<<⎪ ⎪∂∂⎪⎝⎭⎨==-≤≤⎪⎪==≤≤⎩二、用椭圆型第一边值问题的五点差分格式得到线性方程组为2,1,1,,1,10,1,,0,141,?,?,?,?0,1i j i j i j i j i j ijj N j i i N u u u u u h f i j N u u u u i j N -+-+++----=≤≤====≤≤+, 写成矩阵形式Au=f 。

其中1.三 、编写求解线性方程组Au=f 的算法程序, 用下列方法编程计算, 并比较计算速度。

2.用Jacobi 迭代法求解线性方程组Au=f 。

3.用块Jacobi 迭代法求解线性方程组Au=f 。

4. 用SOR 迭代法求解线性方程组Au=f,用试算法确定最佳松弛因子。

1122N N v b v b u f v b ⎛⎫⎛⎫ ⎪ ⎪ ⎪ ⎪== ⎪ ⎪ ⎪ ⎪⎝⎭⎝⎭4114114ii A -⎛⎫ ⎪- ⎪= ⎪- ⎪-⎝⎭11,12,1,121,22,2,21,2,,2211,12,1,121,22,2,221,2,,(,,...,),(,,...,),......,(,,...,)(,,...,)?,(,,...,)?,......,(,,...,)?1,999,0.10.011T T N N TN N N N N T T N N T N N N N N v u u u v u u u v u u u b h f f f b h f f f b h f f f h N h N ====+=+=+===+取或则或,1,,1,2,...,i j f i j N== 1122NN A I I A A I I A -⎛⎫ ⎪- ⎪= ⎪- ⎪-⎝⎭5.用块SOR 迭代法求解线性方程组Au=f,用试算法确定最佳松弛因子。

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

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

数值分析上机实验报告x k x k - f(X k) f (X k)《数值分析》上机实验报告1. 用Newt on法求方程X7-X4+14=0在(0.1,1.9)中的近似根(初始近似值取为区间端点,迭代6次或误差小于0.00001 )。

1.1理论依据:设函数在有限区间[a,b]上二阶导数存在,且满足条件1. f(x)f(b) 02. f(x)在区间[a, b]上不变号3f(x) = 0;4」f (c)〔f .(x) |,其中c是a,b中使mir(| f .(a), f .(b) |)达到的一个b -a则对任意初始近似值x0• [a,b],由Newton迭代过程込f(x k )X“ M(Xk) = Xk — T^,k = 0,1,2,3…f'(X k)所生的迭代序列 % [平方收敛于方程f(x)=0在区间[a,b]上的惟一解: 令7 4f(x)=x -28x 14, f (0.1) 0, f(1.9) ::0f (x) =7x6-112x3=7x3(x3-16) ::: 0f (x) =42x5-336x2=42x2(x3-8) :: 0f (1.9) f (1.9) 0故以1.9为起点x0 =1.9如此一次一次的迭代,逼近X的真实根。

当前后两个的差<=出寸,就认为求出了近似的根。

本程序用Newton法求代数方程(最高次数不大于10)在(a,b )区间的根//限制循环次数1.2 C 语言程序原代码:#i nclude<stdio.h> #in clude<math.h> mai n() {double x2,f,f1; double x1=1.9; // 取初值为 1.9do {x2=x1;f=pow(x2,7)-28*pow(x2,4)+14; f1=7*pow(x2,6)-4*28*pow(x2,3); x 仁 x2-f/f1;}while(fabs(x1-x2)>=0.00001||x1<0.1); printf("计算结果:x=%f\n",x1);}1.3运行结果:* D:\VC + +\EXERCIS E\Debu g\l1.4 MATLAB上机程序fun cti on y=Newt on( f,df,x0,eps,M)d=0;for k=1:Mif feval(df,x0)==0d=2; breakelsex1=x0-feval(f,x0)/feval(df,x0);ende=abs(x1-x0);x0=x1;if e<=eps&&abs(feval(f,x1))v=epsd=1; breakendendif d==1y=x1;elseif d==0y='迭代M次失败';elsey=奇异’endfun cti on y=df(x)y=7*x A6-28*4*x A3;Endfunction y=f(x) y=x A7-28*x A4+14;End>> x0=1.9;>> eps=0.00001;>> M=100;>> x=Newto n('f,'df,x0,eps,M);>> vpa(x,7)1.5问题讨论: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;%%%误差R=1;k=0;%迭代次数初值while (R>5e-6) ;c=(a+b)/2;if f12(a)*f12(c)>0;a=c;elseb=c;endR=b-a;%求出误差k=k+1;endx=c%给出解Newton法及改进的Newton法源程序:clear%%%% 输⼊函数f=input('请输⼊需要求解函数>>','s')%%%求解f(x)的导数df=diff(f);%%%改进常数或重根数miu=2;%%%初始值x0x0=input('input initial value x0>>');k=0;%迭代次数max=100;%最⼤迭代次数R=eval(subs(f,'x0','x'));%求解f(x0),以确定初值x0时否就是解while (abs(R)>1e-8)x1=x0-miu*eval(subs(f,'x0','x'))/eval(subs(df,'x0','x'));R=x1-x0;x0=x1;k=k+1;if (eval(subs(f,'x0','x'))<1e-10);breakendif k>max;%如果迭代次数⼤于给定值,认为迭代不收敛,重新输⼊初值ss=input('maybe result is error,choose a new x0,y/n?>>','s');if strcmp(ss,'y')x0=input('input initial value x0>>');k=0;elsebreakendendendk;%给出迭代次数x=x0;%给出解结果分析和讨论:1.⽤⼆分法计算⽅程在[1,2]内的根。

相关主题
  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

一.实验任务
用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)
%LEASTP.m:least-square fitting with legendre polynomials
%func 指被逼近函数,调用需要用句柄
%a,b 分别指被逼近函数的区间上下限
%n 指最佳平方逼近的次数
syms t;
syms x;
%以Lengendre 多项式为基,构造任意次数的最佳平方逼近多项式
p(2)=t;
p(1)=1;
if n>1
for j=3:1:(n+1)
p(j)=((2*j-3)*t*p(j-1)-(j-2)*p(j-2))/(j-1);
end
end
%变量代换,区间调整为[-1,1]
f=feval(func,(b-a)/2*t+(b+a)/2);
%计算组合系数,其中调用变步长复化梯形求积函数trapz
for j=1:1:(n+1)
c(j)=(2*j-1)/2*trapz(f*p(j),-1,1);
end
%将组合系数与对应的最佳平方多项式相乘然后求和,得到最佳逼近函数sx=0;
for j=1:1:(n+1)
sx=sx+c(j)*p(j);
end
%将变量替换还原
sx=subs(sx,(2*x-a-b)/(b-a));
%使用fplot绘制原函数图像
f1=feval(func,x);
f1=inline(f1);
[x,y]=fplot(f1,[a,b]);
plot(x,y,'r-','linewidth',1.5);
hold on;
%使用plot绘制最佳平方逼近函数图像
g=linspace(a,b,(b-a)*300);
fsx=subs(sx,g);
plot(g,fsx,'b-','linewidth',1.5);
str=strcat(num2str(n),'次最佳平方逼近');
legend('原函数',str);
end
2. 计算组合系数,变步长复化梯形求积法
function To1=trapz(func,a,b)
%半分区间复化梯形公式计算定积分
%func指需要求积分的原函数
%a,b分别指积分上下区间
%初值
h=b-a;
To=(subs(func,a)+subs(func,b))*(b-a)/2;
e=1;
while e>10^-6%迭代终止条件,前后两次积分值差小于10^-6 H=0;
x=a+h/2;
while x<b
H=H+subs(func,x);%计算出所有二分新出现的值的和
x=x+h;
end
To1=0.5*(To+h*H);%计算出新的积分值
e=abs(To1-To);
h=h/2;%继续半分区间,进行迭代计算
To=To1;
end
end
3. 以.m文件定义被逼近函数
function y=func1(x)
y=x*cos(x);
end
五.实验结果
1. 一次最佳平方逼近
c =
-1.1702 -2.4235
sx=
1.253290 - 1.211752*x
2. 二次最佳平方逼近
c =
-1.1702 -2.4235 -0.4265
sx=
-0.159939*x^2 - 0.571997*x + 0.826787
3. 三次最佳平方逼近
c =
-1.1702 -2.4235 -0.4265 1.2216
sx=
0.381759*x^3 - 2.450495*x^2 + 3.092892*x - 0.394843
4. 四次最佳平方逼近
c =
-1.1702 -2.4235 -0.4265 1.2216 0.3123
sx =
0.085392*x^4 - 0.301375*x^3 - 0.693864*x^2 + 1.531443*x - 0.082553
六.分析与讨论
从次数从1到4的最佳平方逼近图像对比可以发现,次数越高,图像拟合效果越好。

尤其是当n=4时,几乎与原函数重合,精度非常高。

但是,从实验中也发现,逼近次数越高,运算量越大,耗时也越长,换言之,时间和空间复杂度越高。

因此,在实际运用中应该从精度要求和算法复杂度两方面来综合考虑最佳平方逼近的次数。

相关文档
最新文档