计算方法上机实验报告
数值分析(计算方法)课程设计实验报告(附程序)
n=4 时,max[L(X)-h(X)]=0.4020;
n=8 时,max[L(X)-h(X)]=0.1708;
n=10 时,max[L(X)-h(X)]=0.1092。
图象分析: 从图象可以看出随着插值节点数的增加出现异常的摆动,中间能较好的接近 原函数,但两边却出现很大的误差。
(3).对定义在(-5,5)上的函数
程序代码 2:
x=[-1:0.2:1]; y=1./(1+25.*x.^2); x0=[-1:0.01:1]; y0=lagrange(x,y,x0); y1=1./(1+25.*x0.^2);
plot(x0,y0,'--r'); hold on; plot(x0,y1,'-b'); x2=abs(y0-y1); max(x2) ; 程序代码3: n=3; for i=1:n x(i)=cos(((2.*i-1).*pi)./(2.*(n+1))); y(i)=1./(1+25.*x(i).*x(i)); end x0=-1:0.01:1; y0=lagrange(x,y,x0); y1=1./(1+25.*x0.^2); plot(x0,y0,'--r') hold on plot(x0,y1,'-b')
以 x1,x2,„,xn+1 为插值节点构造上述各函数的 Lagrange 插值多项式, 比较其 结果。
设计过程: 已知函数 f(x)在 n+1 个点 x0,x1,…,xn 处的函数值为 y0,y1,…,yn 。 求一 n 次多 项式函数 Pn(x),使其满足: Pn(xi)=yi,i=0,1,…,n. 解决此问题的拉格朗日插值多项式公式如下
上机实验报告(精选11篇)
上机实验报告篇1用户名se××××学号姓名学院①实验名称:②实验目的:③算法描述(可用文字描述,也可用流程图):④源代码:(.c的文件)⑤用户屏幕(即程序运行时出现在机器上的画面):2.对c文件的要求:程序应具有以下特点:a可读性:有注释。
b交互性:有输入提示。
c结构化程序设计风格:分层缩进、隔行书写。
3.上交时间:12月26日下午1点-6点,工程设计中心三楼教学组。
请注意:过时不候哟!四、实验报告内容0.顺序表的插入。
1.顺序表的删除。
2.带头结点的单链表的\'插入。
3.带头结点的单链表的删除。
注意:1.每个人只需在实验报告中完成上述4个项目中的一个,具体安排为:将自己的序号对4求余,得到的数即为应完成的项目的序号。
例如:序号为85的同学,85%4=1,即在实验报告中应完成顺序表的删除。
2.实验报告中的源代码应是通过编译链接即可运行的。
3.提交到个人空间中的内容应是上机实验中的全部内容。
上机实验报告篇2一、《软件技术基础》上机实验内容1.顺序表的建立、插入、删除。
2.带头结点的单链表的建立(用尾插法)、插入、删除。
二、提交到个人10m硬盘空间的内容及截止时间1.分别建立二个文件夹,取名为顺序表和单链表。
2.在这二个文件夹中,分别存放上述二个实验的相关文件。
每个文件夹中应有三个文件(.c文件、.obj文件和.exe文件)。
3. 截止时间:12月28日(18周周日)晚上关机时为止,届时服务器将关闭。
三、实验报告要求及上交时间(用a4纸打印)1.格式:《计算机软件技术基础》上机实验报告用户名se××××学号姓名学院①实验名称:②实验目的:③算法描述(可用文字描述,也可用流程图):④源代码:(.c的文件)⑤用户屏幕(即程序运行时出现在机器上的画面):2.对c文件的要求:程序应具有以下特点:a 可读性:有注释。
b 交互性:有输入提示。
数值计算方法上机实验报告
数值计算方法上机实验报告
一、实验目的
本次实验的主要目的是熟悉和掌握数值计算方法,学习梯度下降法的
原理和实际应用,熟悉Python语言的编程基础知识,掌握Python语言的
基本语法。
二、设计思路
本次实验主要使用的python语言,利用python下的numpy,matplotlib这两个工具,来实现数值计算和可视化的任务。
1. 首先了解numpy的基本使用方法,学习numpy的矩阵操作,以及numpy提供的常见算法,如矩阵分解、特征值分解等。
2. 在了解numpy的基本操作后,可以学习matplotlib库中的可视化
技术,掌握如何将生成的数据以图表的形式展示出来。
3. 接下来就是要学习梯度下降法,首先了解梯度下降法的主要原理,以及具体的实际应用,用python实现梯度下降法给出的算法框架,最终
可以达到所期望的优化结果。
三、实验步骤
1. 熟悉Python语言的基本语法。
首先是熟悉Python语言的基本语法,学习如何使用Python实现变量
定义,控制语句,函数定义,类使用,以及面向对象编程的基本概念。
2. 学习numpy库的使用方法。
其次是学习numpy库的使用方法,学习如何使用numpy库构建矩阵,学习numpy库的向量,矩阵操作,以及numpy库提供的常见算法,如矩阵分解,特征值分解等。
3. 学习matplotlib库的使用方法。
计算方法实验报告册
实验一——插值方法实验学时:4实验类型:设计 实验要求:必修一 实验目的通过本次上机实习,能够进一步加深对各种插值算法的理解;学会使用用三种类型的插值函数的数学模型、基本算法,结合相应软件(如VC/VB/Delphi/Matlab/JAVA/Turbo C )编程实现数值方法的求解。
并用该软件的绘图功能来显示插值函数,使其计算结果更加直观和形象化。
二 实验内容通过程序求出插值函数的表达式是比较麻烦的,常用的方法是描出插值曲线上尽量密集的有限个采样点,并用这有限个采样点的连线,即折线,近似插值曲线。
取点越密集,所得折线就越逼近理论上的插值曲线。
本实验中将所取的点的横坐标存放于动态数组[]X n 中,通过插值方法计算得到的对应纵坐标存放于动态数组[]Y n 中。
以Visual C++.Net 2005为例。
本实验将Lagrange 插值、Newton 插值和三次样条插值实现为一个C++类CInterpolation ,并在Button 单击事件中调用该类相应函数,得出插值结果并画出图像。
CInterpolation 类为 class CInterpolation { public :CInterpolation();//构造函数CInterpolation(float *x1, float *y1, int n1);//结点横坐标、纵坐标、下标上限 ~ CInterpolation();//析构函数 ………… …………int n, N;//结点下标上限,采样点下标上限float *x, *y, *X;//分别存放结点横坐标、结点纵坐标、采样点横坐标float *p_H,*p_Alpha,*p_Beta,*p_a,*p_b,*p_c,*p_d,*p_m;//样条插值用到的公有指针,分别存放i h ,i α,i β,i a ,i b ,i c ,i d 和i m};其中,有参数的构造函数为CInterpolation(float *x1, float *y1, int n1) {//动态数组x1,y1中存放结点的横、纵坐标,n1是结点下标上限(即n1+1个结点) n=n1;N=x1[n]-x1[0]; X=new float [N+1]; x=new float [n+1]; y=new float [n+1];for (int i=0;i<=n;i++) {x[i]=x1[i]; y[i]=y1[i]; }for (int i=0;i<=N;i++) X[i]=x[0]+i; }2.1 Lagrange 插值()()nn i i i P x y l x ==∑,其中0,()nj i j j ni jx x l x x x =≠-=-∏对于一个自变量x ,要求插值函数值()n P x ,首先需要计算对应的Lagrange 插值基函数值()i l x float l(float xv,int i) //求插值基函数()i l x 的值 {float t=1;for (int j=0;j<=n;j++) if (j!=i)t=t*(xv-x[j])/(x[i]-x[j]); return t; }调用函数l(float x,int i),可求出()n P xfloat p_l(float x) //求()n P x 在一个点的插值结果 {float t=0;for (int i=0;i<=n;i++) t+=y[i]*l(x,i); return t; }调用p_l(float x)可实现整个区间的插值float *Lagrange() //求整个插值区间上所有采样点的插值结果 {float *Y=new float [N+1]; for (int k=0;k<=N;k++) Y[k]=p_l(x[0]+k*h); return Y; } 2.2Newton 插值010()(,,)()nn i i i P x f x x x x ω==∑,其中101,0()(),0i i j j i x x x i ω-==⎧⎪=⎨-≠⎪⎩∏,0100,()(,,)()ik i nk k j j j kf x f x x x x x ==≠=-∑∏对于一个自变量x ,要求插值函数值()n P x ,首先需要计算出01(,,)i f x x x 和()i x ωfloat *f() {//该函数的返回值是一个长度为n +1的动态数组,存放各阶差商 }float w(float x, int i) {//该函数计算()i x ω }在求()n P x 的函数中调用*f()得到各阶差商,然后在循环中调用w(float x)可得出插值结果 float p_n(float x) {//该函数计算()n P x 在一点的值 }调用p_n(float x)可实现整个区间的插值 float *Newton() {//该函数计算出插值区间内所有点的值 }2.3 三次样条插值三次样条插值程序可分为以下四步编写: (1) 计算结点间的步长i hi 、i α、i β;(2) 利用i hi 、i α、i β产生三对角方程组的系数矩阵和常数向量; (3) 通过求解三对角方程组,得出中间结点的导数i m ; (4) 对自变量x ,在对应区间1[,]i i x x +上,使用Hermite 插值; (5)调用上述函数,实现样条插值。
东南大学计算方法实验报告
计算方法与实习实验报告学院:电气工程学院指导老师:***班级:160093******学号:********实习题一实验1 拉格朗日插值法一、方法原理n次拉格朗日插值多项式为:L n(x)=y0l0(x)+y1l1(x)+y2l2(x)+…+y n l n(x)n=1时,称为线性插值,L1(x)=y0(x-x1)/(x0-x1)+ y1(x-x0)/(x1-x0)=y0+(y1-x0)(x-x0)/(x1-x0)n=2时,称为二次插值或抛物线插值,精度相对高些L2(x)=y0(x-x1)(x-x2)/(x0-x1)/(x0-x2)+y1(x-x0)(x-x2)/(x1-x0)/(x1-x2)+y2(x-x0)(x-x1)/(x2-x0)/(x2-x1)二、主要思路使用线性方程组求系数构造插值公式相对复杂,可改用构造方法来插值。
对节点x i(i=0,1,…,n)中任一点x k(0<=k<=n)作一n 次多项式l k(x k),使它在该点上取值为1,而在其余点x i(i=0,1,…,k-1,k+1,…,n)上为0,则插值多项式为L n(x)=y0l0(x)+y1l1(x)+y2l2(x)+…+y n l n(x) 上式表明:n 个点x i(i=0,1,…,k-1,k+1,…,n)都是l k(x)的零点。
可求得l k三.计算方法及过程:1.输入节点的个数n2.输入各个节点的横纵坐标3.输入插值点4.调用函数,返回z函数语句与形参说明程序源代码如下:#include<iostream>#include<math.h>using namespace std;#define N 100double fun(double *x,double *y, int n,double p);void main(){int i,n;cout<<"输入节点的个数n:";cin>>n;double x[N], y[N],p;cout<<"please input xiangliang x= "<<endl;for(i=0;i<n;i++)cin>>x[i];cout<<"please input xiangliang y= "<<endl;for(i=0;i<n;i++)cin>>y[i];cout<<"please input LagelangrichazhiJieDian p= "<<endl;cin>>p;cout<<"The Answer= "<<fun(x,y,n,p)<<endl;system("pause") ;}double fun(double x[],double y[], int n,double p){double z=0,s=1.0;int k=0,i=0;double L[N];while(k<n){ if(k==0){ for(i=1;i<n;i++)s=s*(p-x[i])/(x[0]-x[i]);L[0]=s*y[0];k=k+1;}else{s=1.0;for(i=0;i<=k-1;i++)s=s*((p-x[i])/(x[k]-x[i]));for(i=k+1;i<n;i++) s=s*((p-x[i])/(x[k]-x[i]));L[k]=s*y[k];k++;}}for(i=0;i<n;i++)z=z+L[i];return z;}五.实验分析n=2时,为一次插值,即线性插值n=3时,为二次插值,即抛物线插值n=1,此时只有一个节点,插值点的值就是该节点的函数值n<1时,结果都是返回0的;这里做了n=0和n=-7两种情况3<n<100时,也都有相应的答案常用的是线性插值和抛物线插值,显然,抛物线精度相对高些n次插值多项式Ln(x)通常是次数为n的多项式,特殊情况可能次数小于n.例如:通过三点的二次插值多项式L2(x),如果三点共线,则y=L2(x)就是一条直线,而不是抛物线,这时L2(x)是一次式。
计算方法与计算 实验一误差分析
% 输出的量--每次迭代次数k和迭代值xk,
%
--每次迭代的绝对误差juecha和相对误差xiangcha,
误差分析
误差问题是数值分析的基础,又是数值分析中一个困难的课题。在实际计算 中,如果选用了不同的算法,由于舍入误差的影响,将会得到截然不同的结果。 因此,选取算法时注重分析舍入误差的影响,在实际计算中是十分重要的。同时, 由于在数值求解过程中用有限的过程代替无限的过程会产生截断误差,因此算法 的好坏会影响到数值结果的精度。 一、实验目的
因为运行后输出结果为: y 1.370 762 168 154 49, yˆ =1.370 744 664 189
38, R 1.750 396 510 491 47e-005, WU= 1.782 679 830 970 664e-005 104 . 所
以, yˆ 的绝对误差为 10 4 ,故 y
③ 运行后输出计算结果列入表 1–1 和表 1-2 中。
④ 将算法 2 的 MATLAB 调用函数程序的函数分别用 y1=15-2*x^2 和
y1=x-(2*x^2+x-15)/(4*x+1)代替,得到算法 1 和算法 3 的调用函数程序,将其保
存,运行后将三种算法的前 8 个迭代值 x1, x2 ,, x8 列在一起(见表 1-1),进行
的精确解 x* 2.5 比较,观察误差的传播.
算法 1 将已知方程化为同解方程 x 15 2x2 .取初值 x0 2 ,按迭代公式
xk1 15 2xk2
C语言上机实验
实验一(第1章实验)实验目的:1.掌握运行C语言程序的全过程。
2.熟悉编译环境。
3.初步熟悉C语言程序的语法规定。
4.了解简单函数的使用方法。
实验内容:1.编程且上机运行:求3个整数的和。
2.编程且上机运行:求2个数的和、差、积和商。
3.编程且上机运行:输入3个数,求最大值。
4.编程且上机运行:输入圆的半径,求圆的面积和周长。
5.在屏幕上输出:“hello world!”实验结果:实验二(第3章实验)1.实验目的:理解C语言的类型系统。
实验内容:写程序测试数据-2在类型char,int,unsigned int,long int,unsigned long int 中存储情况。
实验过程:实验结果:参见各种类型的存储实现描述。
2.实验目的:了解混合类型计算中类型的转换规则。
实验内容:写程序测试多种类型数据一起运算时类型的转换及表达式结果的类型。
注意unsigned int和int数据运算时类型转换的方向。
实验过程:/** 类型转换问题* 试问下面两个表达式等价吗?*/#include <stdio.h>#include <stdlib.h>int main() {unsigned int ui,uj;ui = 1;uj = 2;if (ui < uj)printf("\n%u < %u is true !\n", ui, uj);elseprintf("\n%u < %u is false !\n", ui, uj);if (ui - uj < 0)printf("\n%u - %u <0 is true !\n", ui, uj);elseprintf("\n%u - %u <0 is false !\n", ui, uj);system("pause");return 0;}实验结果:参见类型转换规则。
计算方法二分法实验
计算方法上机实验报告土木0804上机题目;利用二分法求下列方程在【5.5-6.5】的根实验目的:本实验通过二分法对方程求解,通过观察方程x1,x2随二分次数增加而起的变化,以及通过改变方程中七次项系数微小变化对方程求解的影响,加深对二分法求解方程的理解。
并锻炼C语言编程能力。
实验分析:二分法求解采用的是根的存在定理(零点定理):f(x)为[a,b]上的连续函数,若f(a)·f(b)<0,则[a,b]中至少有一个实根。
如果f(x)在[a,b]上还是单调递增或递减的,则f(x)=0仅有一个实根。
取[a,b]的中点x=(a +b)/2作为问题的近似解,那么我们可以估计出绝对误差限仅为区间长的一半,即e=(b-a)/2。
如果这个结果能满足精度要求,我们就停止进一步的计算;如果不能,就求出f(x),结果只能是下面三种情况之一:(1) f(a)·f(x)<0,此时我们有x*∈[a,x];(2) f(x)·f(b)<0,此时我们有x*∈[x,b];(3) f(x)=0,此时x即为问题的精确解。
在前两种情况下,我们可以用x分别替换原问题中的b或a,从而把求解的区间减小了一半。
这样我们又可以取新区间[a,b]的中点。
经过N次迭代后,剩下的区间长为(b- a)/2N。
这也是结果的绝对误差限。
如此继续下去就达到是有根区间逐步缩小的目的,在这一些相互包含的子区间中构造收敛的数列{}kx来逼近根*x。
实验过程:1;C语言程序# include <stdio.h># include <conio.h># include <math.h># define PRECISION 1e12double f(double c){double m,x,m1,m2,m3;x=c;m1=x*x;m2=m1*x;m3=m1*m1;m=m3*m3-36*m3*m2+546*m2*m3-4536*m3*x+22449*m3-67284*m2+118124*m1-109584*x+40320; return(m);}double Root(float a,float b){ int i;i=0;double x1,x2;x1=a;x2=b;do{printf ("\n i=%d,x1=%f,x2=%f ",i,x1,x2);if(f(x1)*f((x1+x2)/2)<0)x2=(x1+x2)/2;elsex1=(x1+x2)/2;i++;}while (i<=100);return (x1);}void main(){float x1,x2,r;printf ("\n input x1,x2: \n ") ;scanf ("%f %f",&x1,&x2);r=Root(x1,x2);printf ("\n Root x=%f",r);}2.结果分析当a7=-36是,迭代到i=21时,出现x1=x2=6.500000,不再发生变化,出现精确解 x=root=x1=6.500000,收敛速度在数列xn 越靠近根时越慢。
计算方法实验报告2
实验报告2:解线性方程组的直接法姓名:杜娟学号:08012324 班级:勘查08-3班一.上机题目用高斯列主元消去法和LU分解法解线性方程组二.目的要求掌握用高斯列主元消去法和LU分解法设计程序,从而实现解线性方程组。
三.方法原理1.如果在一列中选取按模最大的元素,将其调到主干方程位置再做消元,则称为列主元消元法。
调换方程组的次序是为了使运算中做分母量的绝对值尽量地大,减少舍入误差的影响。
2.由高斯消元法得到启发,对消元的过程相当于将分解为一个上三角矩阵和一个下三角矩阵的过程。
如果直接分解得到和,。
这时方程化为,令,由解出;再由,解出。
这就是直接分解法。
四.算法步骤列主元消元法算法1.输入:方程组阶数n,方程组系数矩阵A和常数向量项b。
2.for k=1 to n-1 //选主元的消元过程{//选择{s=|a kk|,m=kfor u=k+1 to nif |a uk|>s then{m=u,s=| a uk|}for v=k to n //交换第k行和第m行{t=a kv; a kv=a mv; a mv=t}t=b k;b k=b m;b m=t}for i=k+1 to n{t=a ik/a kkfor j=k+1 to n{a ij=a ij-t*a kj}b i=b i-t*a kj}}3.for i:=n TO 1 //回代求解4.输出方程组的解 x i, i=1,2,…,n。
如果对于第k步,从k行至n行和从k列至n列中选取按模最大的,对第行和第行交换,对第列和第v列交换,这就是全主元消元法。
在k列和第v列交换时,还要记录下v的序号,以便恢复未知量xk和xv的位置。
LU分解法1计算的第一行元素要计算,则列出式(3.20)等号两边的第1行第1列元素的关系式:故。
一般地,由的第一行元素的关系式得到2计算的第一列元素要计算,则列出式(3.20)等号两边的第2行第1列元素的关系式:故。
计算方法实验报告(附代码)
实验一 牛顿下山法实验说明:求非线性方程组的解是科学计算常遇到的问题,有很多实际背景.各种算法层出不穷,其中迭代是主流算法。
只有建立有效的迭代格式,迭代数列才可以收敛于所求的根。
因此设计算法之前,对于一般迭代进行收敛性的判断是至关重要的。
牛顿法也叫切线法,是迭代算法中典型方法,只要初值选取适当,在单根附近,牛顿法收敛速度很快,初值对于牛顿迭代 至关重要。
当初值选取不当可以采用牛顿下山算法进行纠正。
牛顿下山公式:)()(1k k k k x f x f x x '-=+λ下山因子 ,,,,322121211=λ下山条件|)(||)(|1k k x f x f <+实验代码:#include<iostream> #include<iomanip> #include<cmath>using namespace std;double newton_downhill(double x0,double x1); //牛顿下山法函数,返回下山成功后的修正初值double Y; //定义下山因子Y double k; //k为下山因子Y允许的最小值double dfun(double x){return 3*x*x-1;} //dfun()计算f(x)的导数值double fun1(double x){return x*x*x-x-1;} //fun1()计算f(x)的函数值double fun2(double x) {return x-fun1(x)/dfun(x);} //fun2()计算迭代值int N; //N记录迭代次数double e; //e表示要求的精度int main(){double x0,x1;cout<<"请输入初值x0:";cin>>x0;cout<<"请输入要求的精度:";cin>>e;N=1;if(dfun(x0)==0){cout<<"f'(x0)=0,无法进行牛顿迭代!"<<endl;}x1=fun2(x0);cout<<"x0"<<setw(18)<<"x1"<<setw(18)<<"e"<<setw(25)<<"f(x1)-f(x0)"<<endl;cout<<setiosflags(ios::fixed)<<setprecision(6)<<x0<<" "<<x1<<" "<<fabs(x1-x0)<<" "<<fabs(fun1(x1))-fabs(fun1(x0))<<endl;if(fabs(fun1(x1))>=fabs(fun1(x0))){ //初值不满足要求时,转入牛顿下山法x1=newton_downhill(x0,x1);} //牛顿下山法结束后,转入牛顿迭代法进行计算while(fabs(x1-x0)>=e){ //当精度不满足要求时N=N+1;x0=x1;if(dfun(x0)==0){cout<<"迭代途中f'(x0)=0,无法进行牛顿迭代!"<<endl;} x1=fun2(x0);cout<<setiosflags(ios::fixed)<<setprecision(6)<<x0<<" "<<x1<<" "<<fabs(x1-x0)<<endl;}cout<<"迭代值为:"<<setiosflags(ios::fixed)<<setprecision(6)<<x1<<'\n';cout<<"迭代次数为:"<<N<<endl;return 0;}double newton_downhill(double x0,double x1){Y=1;cout<<"转入牛顿下山法,请输入下山因子允许的最小值:";cin>>k;while(fabs(fun1(x1))>=fabs(fun1(x0))){if(Y>k){Y=Y/2;}else {cout<<"下山失败!";exit(0);}x1=x0-Y*fun1(x0)/dfun(x0);}//下山成功则cout<<"下山成功!Y="<<Y<<",转入牛顿迭代法计算!"<<endl;return x1;}实验结果:图4.1G-S 迭代算法流程图实验二 高斯-塞德尔迭代法实验说明:线性方程组大致分迭代法和直接法。
数值计算方法上机实验报告
数值计算方法上机实验报告实验目的:复习和巩固数值计算方法的基本数学模型,全面掌握运用计算机进行数值计算的具体过程及相关问题。
利用计算机语言独立编写、调试数值计算方法程序,培养学生利用计算机和所学理论知识分析解决实际问题的能力。
上机练习任务:利用计算机基本C 语言编写并调试一系列数值方法计算通用程序,并能正确计算给定题目,掌握调试技能。
掌握文件使用编程技能,如文件的各类操作,数据格式设计、通用程序运行过程中文件输入输出运行方式设计等。
一、各算法的算法原理及计算机程序框图1. 列主元高斯消去法算法原理:高斯消去法是利用现行方程组初等变换中的一种变换,即用一个不为零的数乘一个方程后加只另一个方程,使方程组变成同解的上三角方程组,然后再自下而上对上三角方程组求解。
列选住院是当高斯消元到第k 步时,从k 列的kk a 以下(包括kk a )的各元素中选出绝对值最大的,然后通过行交换将其交换到kk a 的位置上。
交换系数矩阵中的两行(包括常数项),只相当于两个方程的位置交换了,因此,列选主元不影响求解的结果。
●源程序:#define N 200#include "stdio.h"#include "math.h"FILE *fp1,*fp2;void LZ(){int n,i,j,k=0,l;double d,t,t1;static double x[N],a[N][N];fp1=fopen("a1.txt","r");fp2=fopen("b1.txt","w");fscanf(fp1,"%d",&n);for(i=0;i<n;++i)for(j=0;j<=n;++j){fscanf(fp1,"%lf",&a[i][j]);}{d=a[k][k];l=k;i=k+1;do{if(fabs(a[i][k])>fabs(d)) /*选主元*/{d=a[i][k];l=i;}i++;}while(i<n);if(d==0){printf("\n输入矩阵有误!\n");}else{ /*换行*/if(l!=k){for(j=k;j<=n;j++){t=a[l][j];a[l][j]=a[k][j];a[k][j]=t;}}}for(j=k+1;j<=n;j++) /*正消*/ a[k][j]/=a[k][k];for(i=k+1;i<n;i++)for(j=k+1;j<=n;j++)a[i][j]-=a[i][k]*a[k][j];k++;}while(k<n);if(k!=0){for(i=n-1;i>=0;i--) /*回代*/ {t1=0;for(j=i+1;j<n;j++)t1+=a[i][j]*x[j];x[i]=a[i][n]-t1;}for(i=0;i<n;i++)fprintf(fp2,"\n 方程组的根为x[%d]=%lf",i+1,x[i]); fclose(fp1); fclose(fp2); }main() { LZ(); }● 具体算例及求解结果:用列选主元法求解下列线性方程组⎪⎩⎪⎨⎧=++=++=-+28x x 23x 2232832321321321x x x x x x 输入3 输出结果:方程组的根为x[1]=6.0000001 2 -3 8 方程组的根为x[2]=4.000000 2 1 3 22 方程组的根为x[3]=2.000000 3 2 1 28● 输入变量、输出变量说明:输入变量:ij a 系数矩阵元素,i b 常向量元素 输出变量:12,,n b b b 解向量元素2. 杜里特尔分解法解线性方程● 算法原理:求解线性方程组Ax b =时,当对A 进行杜里特尔分解,则等价于求解LUx b =,这时可归结为利用递推计算相继求解两个三角形(系数矩阵为三角矩阵)方程组,用顺代,由Ly b =求出y ,再利用回带,由Ux y =求出x 。
上机实验报告怎么写
上机实验报告怎么写实验目的本次实验的目的是...实验环境本次实验使用的工具和软件环境如下:- 操作系统:Windows 10- 开发工具:Visual Studio Code- 编程语言:Python实验步骤步骤一:准备工作在开始实验之前,我们需要进行一些准备工作,包括安装相应的软件环境和准备实验材料等。
具体的准备工作如下:1. 安装操作系统:确保使用最新的操作系统,并安装所需的驱动程序。
2. 安装开发工具:下载并安装Visual Studio Code,并配置相应的插件和设置。
3. 安装Python:下载并安装Python解释器,并配置环境变量。
步骤二:实验设计在本次实验中,我们将实现一个简单的计算器功能,包括加法、减法、乘法和除法。
具体的实验设计如下:1. 设计界面:使用Tkinter库创建一个窗口,并在窗口中添加输入框和按钮等组件。
2. 实现功能:根据用户的输入,进行相应的数值计算,并将结果显示在窗口中。
步骤三:编写代码根据实验设计,我们开始编写代码。
代码的具体实现如下:pythonimport tkinter as tkdef calculate():获取用户输入的数值num1 = float(entry1.get())num2 = float(entry2.get())执行计算操作result = num1 + num2显示计算结果label.config(text="计算结果:" + str(result))创建窗口window = ()添加输入框entry1 = tk.Entry(window)entry1.pack()entry2 = tk.Entry(window)entry2.pack()添加按钮button = tk.Button(window, text="计算", command=calculate) button.pack()添加结果显示标签label = bel(window, text="计算结果:")label.pack()启动事件循环window.mainloop()步骤四:实验结果在编写完代码后,我们进行了实验测试,并记录下了实验结果。
《计算方法》实验报告材料
double Newton(double x,vector<double>&X,vector<double>&Y);
int main(){
char a='n';
do{
int n;
cout<<"请输入插值点个数:"<<endl;
for(int i=0;i<N;i++){
X[i]=p;
Y[i]=1/(1+p*p);
p=p+c;
}
cout<<"请输入要求值x的值:"<<endl;
double x;
cin>>x;
double result=fenduan(N,X,Y,x,c);
cout<<"由分段线性插值法得出结果: "<<result<<endl;
cin>>n;
vector<double>X(n,0);
vector<double>Y(n,0);
cout<<"请输入插值点对应的值及函数值(Xi,Yi):"<<endl;
for(int i=0;i<n;i++){
cin>>X[i]>>Y[i];
}
cout<<"请输入要求值x的值:"<<endl;
数值分析上机实验报告
一、实验目的通过本次上机实验,掌握数值分析中常用的算法,如二分法、牛顿法、不动点迭代法、弦截法等,并能够运用这些算法解决实际问题。
同时,提高编程能力,加深对数值分析理论知识的理解。
二、实验环境1. 操作系统:Windows 102. 编程语言:MATLAB3. 实验工具:MATLAB数值分析工具箱三、实验内容1. 二分法求方程根二分法是一种常用的求方程根的方法,适用于连续函数。
其基本思想是:从区间[a, b]中选取中点c,判断f(c)的符号,若f(c)与f(a)同号,则新的区间为[a, c],否则为[c, b]。
重复此过程,直至满足精度要求。
2. 牛顿法求方程根牛顿法是一种迭代法,适用于可导函数。
其基本思想是:利用函数在某点的导数值,求出函数在该点的切线方程,切线与x轴的交点即为方程的近似根。
3. 不动点迭代法求方程根不动点迭代法是一种迭代法,适用于具有不动点的函数。
其基本思想是:从初始值x0开始,不断迭代函数g(x)的值,直至满足精度要求。
4. 弦截法求方程根弦截法是一种线性近似方法,适用于可导函数。
其基本思想是:利用两点间的直线近似代替曲线,求出直线与x轴的交点作为方程的近似根。
四、实验步骤1. 二分法求方程根(1)编写二分法函数:function [root, error] = bisection(a, b, tol)(2)输入初始区间[a, b]和精度要求tol(3)调用函数计算根:[root, error] = bisection(a, b, tol)2. 牛顿法求方程根(1)编写牛顿法函数:function [root, error] = newton(f, df, x0, tol)(2)输入函数f、导数df、初始值x0和精度要求tol(3)调用函数计算根:[root, error] = newton(f, df, x0, tol)3. 不动点迭代法求方程根(1)编写不动点迭代法函数:function [root, error] = fixed_point(g, x0, tol)(2)输入函数g、初始值x0和精度要求tol(3)调用函数计算根:[root, error] = fixed_point(g, x0, tol)4. 弦截法求方程根(1)编写弦截法函数:function [root, error] = secant(f, x0, x1, tol)(2)输入函数f、初始值x0和x1,以及精度要求tol(3)调用函数计算根:[root, error] = secant(f, x0, x1, tol)五、实验结果与分析1. 二分法求方程根以方程f(x) = x^2 - 2 = 0为例,输入初始区间[a, b]为[1, 3],精度要求tol 为1e-6。
计算方法与实习的实验报告范文
计算方法与实习的实验报告范文1舍入误差与数值稳定性1.1目的与要求(1)通过上机编程,复习巩固以前所学程序设计语言及上机操作指令;(2)通过上机计算,了解舍入误差所引起的数值不稳定性。
1.2舍入误差和数值稳定性1.2.1概要舍入误差在计算方法中是一个很重要的概念。
在实际计算中如果选用了不同的算法,由于舍入误差的影响,将会得到截然不同的结果。
因此,选取稳定的算法在实际计算中是十分重要的。
1.2.2程序和实例对n=0,1,2,…,40计算定积分1某n某5n0d某。
算法利用递推公式yn=ln6-ln50.182322。
程序如下:#include#includevoidmain(){doubley_0=log(6.0/5.0),y_1;intn=1;printf(\while(1){y_1=1.0/n-5某y_0;printf(\if(n>=40)break;y_0=y_1;n++;if(n%2==0)printf(\}}1n5y(n=1,2,…,40)取y0=11某50d某=2方程求根2.1实验目的(1)通过对二分法与牛顿迭代法作编程练习与上级运算,进一步体会二分法与牛顿迭代法的不同特点;(2)编写割线迭代法的程序,求非线性迭代法的解,并与牛顿迭代法作比较。
2.2二分法2.2.1算法给定区间[a,b],并设f(a)与f(b)符号相反,取ε为根的容许误差,δ为|f(某)|的容许误差。
1令c=(a+b)/2;2如果(c-a)0,则令a=c;否则令b=c,重复1,2,3。
2.2.2程序与实例求方程f(某)=某4某100在1.5附近的根。
程序如下:#include#include#defineep5e-6#definedelta1e-6floatBiection(floata,floatb,float(某f)(float)){floatc,fc,fa=(某f)(a),fb=(某f)(b);intn=1;printf(\二分次数\\t\\tc\\t\\tf(c)\\n\while(1){if(fa某fb>0){printf(\不能用二分法求解\c=(a+b)/2,fc=(某f)(c);printf(\if(fab(fc)returnc;}floatf(float某){return某某某某某+4某某某某-10;32}intmain(){floata=1,b=2;float某;某=Biection(a,b,f);printf(\方程的根为%f\return0;}3线性方程组数值解法3.1目的与要求(1)熟悉求解线性方程组的有关理论和方法;(2)会编制列主元消去法、LU分解法、雅克比及高斯赛德尔迭代法的程序;(3)通过实际计算,进一步了解各种方法的优缺点,选择合适的数值方法。
计算方法实验上机报告(完整版)
简单迭代法#include<stdio.h>#include<math.h>#define x0 3.0#define MAXREPT 1000#define EPS 1E-6#define G(x) pow(12*x+sin(x)-1,1.0/3)void main(){int i;double x_k=x0,x_k1=x0;printf("k\txk\n");for(i=0;i<MAXREPT;i++){printf("%d\t%g\n",i,x_k1);x_k1=G(x_k);if (fabs(x_k1-x_k)<EPS){printf("THE ROOT IS x=%g,k=%d\n",x_k1,i);return;}x_k=x_k1;}printf("AFTER %d repeate,no solved.\n",MAXREPT);}结果牛顿迭代法一#include<stdio.h>#include<math.h>#define x0 3.0#define MAXREPT 1000#define EPS 1E-6#define G(x) x-(pow(x,3)-sin(x)-12*x+1)/(3*pow(x,2)-cos(x)-12) void main(){int i;double x_k=x0,x_k1=x0;printf("k\txk\n");for(i=0;i<MAXREPT;i++){printf("%d\t%g\n",i,x_k1);x_k1=G(x_k);if (fabs(x_k1-x_k)<EPS){printf("THE ROOT IS x=%g,k=%d\n",x_k1,i);return;}x_k=x_k1;}printf("AFTER %d repeate,no solved.\n",MAXREPT);}结果埃特金加速法#include<stdio.h>#include<math.h>#define x0 3.0#define MAXREPT 1000#define EPS 1E-6#define G(x) (pow(x,3)-sin(x)+1)/12void main(){int i;double x1=x0,x2=x0;double z,y;printf("k\tx1\tx2\txk\n");for(i=0;i<MAXREPT;i++){if(i==0)printf("%d\t\t\t%g\n",i,x2);elseprintf("%d\t%g\t%g\t%g\n",i,y,z,x2);y=G(x1);z=G(y);x2=z-((z-y)*(z-y))/(z-2*y+x1);if (fabs(x2-x1)<EPS){printf("THE ROOT IS x=%g,k=%d\n",x2,i);return;}x1=x2;}printf("AFTER %d repeate,no solved.\n",MAXREPT);} 结果牛顿迭代法二#include<stdio.h>#include<math.h>#define x0 1.5#define MAXREPT 1000#define EPS 1E-6#define G(x) x-(pow(x,3)+pow(x,2)-3*x-3)/(3*pow(x,2)+2*x-3) void main(){int i;double x_k=x0,x_k1=x0;printf("k\txk\n");for(i=0;i<MAXREPT;i++){printf("%d\t%g\n",i,x_k1);x_k1=G(x_k);if (fabs(x_k1-x_k)<EPS){printf("THE ROOT IS x=%g,k=%d\n",x_k1,i);return;}x_k=x_k1;}printf("AFTER %d repeate,no solved.\n",MAXREPT);}结果弦截法#include<stdio.h>#include<math.h>#define x0 0#define x1 1.5#define MAXREPT 1000#define EPS 1E-6#define G(x) pow(x,3)+pow(x,2)-3*x-3void main(){int i;double x_k=x0,x_k1=x1,x_k2=0;double y,z;printf("k\txk\n");for(i=0;i<MAXREPT;i++){printf("%d\t%g\n",i,x_k2);y=G(x_k);z=G(x_k1);x_k2=x_k1-(z*(x_k1-x_k))/(z-y);if (fabs(x_k2-x_k1)<EPS){printf("THE ROOT IS x=%g,k=%d\n",x_k2,i);return;}x_k=x_k1;x_k1=x_k2;}printf("AFTER %d repeate,no solved.\n",MAXREPT); } 结果高斯顺序消元法#include<stdio.h>#include<math.h>#define N 4static double aa[N][N+1]={{2,4,0,1,1},{3,8,2,2,3},{1,3,3,0,6},{2,5,2,2,3}}; int gauss(double a[][N+2],double x[]);void putout(double a[][N+2]);void main(){int i,j,det;double a[N+1][N+2],x[N+1];for(i=1;i<=N;i++)for(j=1;j<=N+1;j++)a[i][j]=aa[i-1][j-1];det=gauss(a,x);if(det!=0)for(i=1;i<=N;i++)printf(" x[%d]=%g",i,x[i]);printf("\n");}int gauss(double a[][N+2],double x[]){int i,j,k;double c;putout(a);for(k=1;k<=N-1;k++){ if(fabs(a[k][k])<1e-17){printf("\n pivot element is 0.fail!\n");return 0;}for(i=k+1;i<=N;i++){c=a[i][k]/a[k][k];for(j=k;j<=N+1;j++){a[i][j]=a[i][j]-c*a[k][j];}}putout(a);}if(fabs(a[N][N])<1e-17){printf("\n pivot element is 0.fail!\n");return 0;}for(k=N;k>=1;k--){x[k]=a[k][N+1];for(j=k+1;j<=N;j++){x[k]=x[k]-a[k][j]*x[j];}x[k]=x[k]/a[k][k];}return(1);}void putout(double a[][N+2]){for(int i=1;i<=N;i++){for(int j=1;j<=N+1;j++)printf("%-15g",a[i][j]);printf("\n");}printf("\n");}结果雅克比迭代法#include<stdio.h>#include<math.h>#define N 5#define EPS 0.5e-4static double aa[N][N]={{4,-1,0,-1,0},{-1,4,-1,0,-1},{0,-1,4,-1,0},{-1,0,-1,4,-1},{0,-1,0,-1,4}}; static double bb[N]={2,1,2,1,2};void main(){int i,j,k,NO;double a[N+1][N+1],b[N+1],x[N+1],y[N+1];double d,sum,max;for(i=1;i<=N;i++){for(j=1;j<=N;j++)a[i][j]=aa[i-1][j-1];b[i]=bb[i-1];}printf("\n 请输入最大迭代次数(尽量取大值!)NO:");scanf("%d",&NO);printf("\n");for(i=1;i<=N;i++)x[i]=0;k=0;printf(" k",' ');for(i=1;i<=N;i++)printf("%8cx[%d]",' ',i);printf("\n 0");for(i=1;i<=N;i++)printf("%12.8g",x[i]);printf("\n");do{for(i=1;i<=N;i++){sum=0.0;for(j=1;j<=N;j++)if(j!=i) sum=sum+a[i][j]*x[j];y[i]=(-sum+b[i])/a[i][i];}max=0.0;for(i=0;i<=N;i++){d=fabs(y[i]-x[i]);if(max<d) max=d;x[i]=y[i];}printf("%6d",k+1);for(i=1;i<=N;i++)printf("%12.8g",x[i]);printf("\n");k++;}while((max>=EPS)&&(k<NO));printf("\nk=%d\n",k);if(k>=NO) printf("\nfail!\n");elsefor(i=1;i<=N;i++)printf("x[%d]=%g\t",i,x[i]);}结果拉格朗日插值多项式#include<stdio.h>#include<math.h>#define N 4doublex[N]={0.56160,0.56280,0.56401,0.56521},y[N]={0.82741,0.82659,0.82577,0.82495}; void main(){double x=0.5635;double L(double xx);double lagBasis(int k,double xx);void output();output();printf("\n近似值L(%g)=%g\n",x,L(x));}double lagBasis(int k,double xx){double lb=1;int i;for(i=0;i<N;i++)if(i!=k) lb*=(xx-x[i])/(x[k]-x[i]);return lb;}double L(double xx){double s=0;int i;for(i=0;i<=N;i++)s+=lagBasis(i,xx)*y[i];return s;}void output(){int i;printf("\n各节点信息:\nxi:");for(i=0;i<N;i++)printf("\t%g",x[i]);printf("\nyi:");for(i=0;i<N;i++)printf("\t%g",y[i]);}结果牛顿插值多项式#include <math.h>#include <stdio.h>int a;#define M 4double x[M+1]={0.4,0.55,0.65,0.8,0.9},y[M+1]={0.41075,0.57815,0.69675,0.88811,1.02652}; void main(){double x;printf("输入x=");scanf("%lf",&x);printf("次数:");scanf("%d",&a);double N(double xx,int a);void output();output();printf("\n%d次牛顿插值多项式N(%g)=%g\n",a,x,N(x,a));}double N(double xx,int a){double s=y[0],d=1;int i,j;double df[M+1][M+1];for(i=0;i<=M;i++)df[i][0]=y[i];for(j=1;j<=a;j++)for(i=j;i<=a;i++)df[i][j]=(df[i][j-1]-df[i-1][j-1])/(x[i]-x[i-j]);printf("\nx\tf(x)\t");for(j=1;j<=a;j++) printf("%5d阶",j);for(i=0;i<=a;i++){printf("\n%g\t%g",x[i],y[i]);for(j=1;j<=i;j++)printf("\t%7.5g",df[i][j]);}for(i=1;i<=a;i++){d*=(xx-x[i-1]);s+=df[i][i]*d;}return s;}void output(){int i;printf("\n各节点信息:\nxi:");for(i=0;i<=M;i++)printf("\t%7g",x[i]);printf("\nyi:");for(i=0;i<=M;i++)printf("\t%7g",y[i]);}结果复合梯形公式#include<stdio.h>#include<math.h>#define f(x) 1/(x*x+1)#define Pi 3.1415926void main(){double a=0,b=1;double T,h,x;int n,i;printf("please input n:");scanf("%d",&n);h=(b-a)/n;x=a;T=0;for(i=1;i<n;i++){x+=h;T+=f(x);}T=(f(a)+2*T+f(b))*h/2;printf("T(%d)=%g\n",n,T);printf("The exact value is %g\n",Pi/4);}复合辛普森公式#include<stdio.h>#include<math.h>#define f(x) 1/(1+x*x)#define Pi 3.1415926void main(){double a=0,b=1;double S,h,x;int n,i;printf("please input Even n:");scanf("%d",&n);h=(b-a)/n;x=a; S=0;for(i=1;i<n;i++){x+=h;if(i%2==0) S+=2*f(x);else S+=4*f(x);}S=(f(a)+S+f(b))*h/3;printf("S(%d)=%g\n",n,S);printf("The exact value is %g\n",Pi/4);}龙贝格公式加速#include<stdio.h>#include<math.h>#define f(x) sin(x)/(1+x)#define M 3void main(){double a=0,b=1;double Long(double a,double b);printf("近似值I=%g\n",Long(a,b));}double Long(double a,double b){int n=1,i=1,j=1;double T[M+1][M+1],h,x,sum;h=b-a;T[0][0]=(f(a)+f(b))/2;for(j=1;j<=3;j++){x=a;h/=2;n*=2;sum=0;for(i=1;i<=n;i+=2){x=a+i*h;sum+=f(x);}T[j][0]=T[j-1][0]/2+h*sum;}for(i=1;i<=M;i++)for(j=1;j<=i;j++){T[i][j]=(pow(4,j)*T[i][j-1]-T[i-1][j-1])/(pow(4,j)-1);}printf("k\tT0\tT1\tT2\tT3\n");for(i=0;i<=M;i++){printf("%d",i);for(j=0;j<=i;j++)printf(" %g",T[i][j]);printf("\n");}return T[M][M];}。
计算方法实验报告
2019年计算方法(B)实验报告姓名:学号:专业:课程:计算方法(B)目录一、实验综述 (1)二、实验内容 (1)2.1 实验一 (1)2.2 实验二 (2)2.3 实验三 (3)2.4 实验四 (4)2.5 实验五 (6)三、思考总结 (7)附件A1 (8)附件A2 (9)附件A3 (10)附件A4 (12)附件A5 (14)一、实验综述计算方法在工程实践中得到了广泛的应用,是理工类研究生必备的知识技能。
按照2019年计算方法课程学习要求,本文对计算方法上机题目进行了算法设计、分析,利用matlab 2019b版本对算法进行实现,最终形成了实验报告。
以下为本次实验报告具体内容,包括五个实验部分和一个思考总结部分。
二、实验内容2.1 实验一2.1.1 实验题目用Jacobi迭代和Gauss-Seidel迭代解电流方程组,使各部分电流的误差均小于10-3。
2.1.2 算法分析a)首先列出方程组的系数矩阵A以及等式右端的矩阵b,A=[28,-3,0,0,0;-3,38,-10,0,-5;0,-10,25,-15,0;0,0,-15,45,0;0,-5,0,0,30 ];b=[10;0;0;0;0];为了验证A是否收敛,我们通过判断系数矩阵A是否为严格对角占优矩阵进行确定。
如果是,则可以进行Jacobi迭代和Gauss-Seidel迭代(利用matlab程序验证后,证明了矩阵A为严格对角占优矩阵);如果不是,则需要采用其他方法进行判断迭代是否收敛。
b)对矩阵A分裂成三部分,,其中D为A的对角矩阵,E为A的下三角矩阵的相反数,F为A的上三角矩阵的相反数。
c) Jacobi迭代。
取x得初始向量为x=[0;0;0;0;0],利用迭代公式进行循环计算,当的无穷范数小于10-3,即,停止循环。
d) Gauss-Seidel迭代。
取x得初始向量为x=[0;0;0;0;0],利用迭代公式进行循环计算,当的无穷范数小于10-3,即,停止循环。
数值分析上机实验报告
数值分析上机实验报告数值分析上机实验报告一、引言数值分析是一门研究利用计算机进行数值计算的学科。
通过数值分析,我们可以使用数学方法和算法来解决实际问题,例如求解方程、插值和逼近、数值积分等。
本次上机实验旨在通过编程实现数值计算方法,并应用于实际问题中。
二、实验目的本次实验的目的是掌握数值计算方法的基本原理和实现过程,加深对数值分析理论的理解,并通过实际应用提高编程能力。
三、实验内容1. 数值求解方程首先,我们使用二分法和牛顿迭代法分别求解非线性方程的根。
通过编写程序,输入方程的初始值和精度要求,计算得到方程的根,并与理论解进行对比。
2. 数值插值和逼近接下来,我们使用拉格朗日插值和最小二乘法进行数据的插值和逼近。
通过编写程序,输入给定的数据点,计算得到插值多项式和逼近多项式,并绘制出插值曲线和逼近曲线。
3. 数值积分然后,我们使用梯形法和辛普森法进行定积分的数值计算。
通过编写程序,输入被积函数和积分区间,计算得到定积分的近似值,并与解析解进行比较。
四、实验步骤1. 数值求解方程(1)使用二分法求解非线性方程的根。
根据二分法的原理,编写程序实现二分法求解方程的根。
(2)使用牛顿迭代法求解非线性方程的根。
根据牛顿迭代法的原理,编写程序实现牛顿迭代法求解方程的根。
2. 数值插值和逼近(1)使用拉格朗日插值法进行数据的插值。
根据拉格朗日插值法的原理,编写程序实现数据的插值。
(2)使用最小二乘法进行数据的逼近。
根据最小二乘法的原理,编写程序实现数据的逼近。
3. 数值积分(1)使用梯形法进行定积分的数值计算。
根据梯形法的原理,编写程序实现定积分的数值计算。
(2)使用辛普森法进行定积分的数值计算。
根据辛普森法的原理,编写程序实现定积分的数值计算。
五、实验结果与分析1. 数值求解方程通过二分法和牛顿迭代法,我们成功求解了给定非线性方程的根,并与理论解进行了对比。
结果表明,二分法和牛顿迭代法都能够较好地求解非线性方程的根,但在不同的问题中,二者的收敛速度和精度可能会有所差异。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
《计算方法》上机实验报告班级:XXXXXX小组成员:XXXXXXXXXXXXXXXXXXXXXXXXXXXX任课教师:XXX二〇一八年五月二十五日前言通过进行多次的上机实验,我们结合课本上的内容以及老师对我们的指导,能够较为熟练地掌握Newton 迭代法、Jacobi 迭代法、Gauss-Seidel 迭代法、Newton 插值法、Lagrange 插值法和Gauss 求积公式等六种算法的原理和使用方法,并参考课本例题进行了MATLAB 程序的编写。
以下为本次上机实验报告,按照实验内容共分为六部分。
实验一:一、实验名称及题目: Newton 迭代法例2.7(P38):应用Newton 迭代法求x 3−x −1=0在x =1附近的数值解x k ,并使其满足|x k −x k−1|<10−8. 二、解题思路:设'x 是0)(=x f 的根,选取0x 作为'x 初始近似值,过点())(,00x f x 做曲线)(x f y =的切线L ,L 的方程为))((')(000x x x f x f y -+=,求出L 与x 轴交点的横坐标)(')(0001x f x f x x -=,称1x 为'x 的一次近似值,过点))(,(11x f x 做曲线)(x f y =的切线,求该切线与x 轴的横坐标)(')(1112x f x f x x -=称2x 为'x的二次近似值,重复以上过程,得'x 的近似值序列{}n x ,把)(')(1n n n n x f x f x x -=+称为'x 的1+n 次近似值,这种求解方法就是牛顿迭代法。
三、Matlab 程序代码:function newton_iteration(x0,tol) syms z %定义自变量 format long %定义精度 f=z*z*z-z-1;f1=diff(f);%求导 y=subs(f,z,x0);y1=subs(f1,z,x0);%向函数中代值 x1=x0-y/y1; k=1;while abs(x1-x0)>=tol x0=x1;y=subs(f,z,x0); y1=subs(f1,z,x0); x1=x0-y/y1;k=k+1; endx=double(x1) K四、运行结果:实验二:一、实验名称及题目:Jacobi 迭代法例3.7(P74):试利用Jacobi 迭代公式求解方程组[5−1−1−1−110−1−1−1−15−1−1−1−110][x 1x 2x 3x 4]=[−412834]要求数值解X (k)满足||X −X (k )||2≤10−4,其中X =(1,2,3,4)T 为方程组的精确解. 二、解题思路:首先将方程组中的系数矩阵A 分解成三部分,即:U D L A ++=,D 为对角阵,L 为下三角矩阵,U 为上三角矩阵。
之后确定迭代格式,f X B X k k +=+)()1(*,( ⋅⋅⋅=2,1,0k , k 即迭代次数),B 称为迭代矩阵。
最后选取初始迭代向量)0(X ,开始逐次迭代。
最后验证精度。
(迭代阵:b D UXD Xk k1)(1)1(--++-=。
)雅克比迭代法的优点明显,计算公式简单,每迭代一次只需计算一次矩阵和向量的乘法,且计算过程中原始矩阵A 始终不变,比较容易并行计算。
然而这种迭代方式收敛速度较慢,而且占据的存储空间较大。
三、Matlab 程序代码:function jacobi(A,b,x0,eps,x1) D = diag(diag(A));%求A 的对角矩阵 L = -tril(A,-1);%求A 的下三角矩阵U = -triu(A,1);%求A的上三角矩阵B = D\(L+U);f = D\b;x = B*x0+f;n = 1;%迭代次数while norm(x-x1)>=epsx = B*x+f;n = n+1;endformat longnxjingdu=norm(x-x1)四、运行结果:实验三:一、实验名称及题目:Gauss-Seidel 迭代法例 3.8(P75):试利用Gauss-Seidel迭代公式求解方程组[5−1−1−1−110−1−1−1−15−1−1−1−110][x1x2x3x4]=[−412834],并使其数值解X(k)满足精度要求||X−X(k)||2≤10−4,其中X=(1,2,3,4)T为方程组的精确解.二、解题思路:Gauss-Seidel 迭代法与Jacobi 迭代法思路相近,首先将方程组中的系数矩阵A 分解成三部分,即:U D L A ++=,D 为对角阵,L 为下三角矩阵,U 为上三角矩阵。
之后确定迭代格式,f X B X k k +=+)()1(*,( ⋅⋅⋅=2,1,0k , k 即迭代次数),B 称为迭代矩阵。
最后选取初始迭代向量0X ,开始逐次迭代。
最后验证精度。
(迭代阵:b L D UXL D Xk k1)(1)1()()(--++++-=。
)Gauss-Seidel 迭代法与Jacobi 迭代法相比速度更快,但不全如此。
有例子表明:Gauss-Seidel 迭代法收敛时,Jacobi 迭代法可能不收敛;而Jacobi 迭代法收敛时,Gauss-Seidel 迭代法也可能不收敛。
三、Matlab 程序代码:function gauss_seidel(A,b,x0,eps,x1) D = diag(diag(A));%求A 的对角矩阵 L = -tril(A,-1);%求A 的下三角矩阵 U = -triu(A,1);%求A 的上三角矩阵 B = (D-L)\U; f = (D-L)\b; x = B*x0+f;n = 1;%迭代次数while norm(x1-x)>=eps x = B*x+f; n = n+1; endformat long n xjingdu=norm(x1-x)四、运行结果:实验四:一、实验名称及题目: Lagrange 插值法例4.1(P88):给定函数f (x )=x(1+cosx)及插值节点x 0=0,x 1+π8,x 2=π4,x 3=3π8,x 4=π2.试构造Lagrange 插值多项式,给出其误差估计,并由此计算f (3π16)及其误差.二、解题思路:一般来说,如果我们有n 个点()()n n y x y x ,,...,1,1,各i x 互不相同。
那么应用拉格朗日插值公式所得到的拉格朗日插值多项式为:∑==nj j j x l y x L 0)()(,其中每个)(x l j 为拉格朗日基本多项式(或称插值基函数),其表达式为:∏≠=--=nj i i ij ij x x x x x l ,0)(。
三、Matlab程序代码:function y=lagrange(x0,x)n=length(x0);%向量长度s=0;for k=1:n%k从1到n的循环p=1.0;for j=1:nif j~=k%“~=”不等于的意思p=p*(x-x0(j))/(x0(k)-x0(j));endendy0=x0(k)*(1+cos(x0(k)));s=p*y0+s;endformat longswucha=abs(x*(1+cos(x))-s)四、运行结果:五、Lagrange插值图像绘制%Lagrange插值图像算法x=linspace(0,1002,200);s=linspace(0,1000,200);x0=[0;pi/8;pi/4;3*pi/8;pi/2];n=length(x0);s=0;for k=1:np=1.0;for j=1:nif j~=kp=p.*(x-x0(j))/(x0(k)-x0(j));endendy0=x0(k)*(1+cos(x0(k)));s=p*y0+s;endplot(x,s,'r');grid on;title('Lagrange²åֵͼÏñ')xlabel('X'),ylabel('Y');axis normal;实验五:一、实验名称及题目:Newton 插值法例 4.3(P96):已知f (x )=√1+cℎ2x ,试取插值节点x 0=0.35,x 1=0.50,x 2=0.65,x 3=0.80,x 4=0.95,构造4次Newton 插值多项式,由此计算f (0.7)的逼近值,并指出其绝对误差. 二、解题思路:将拉格朗日插值公式∑∑∏==≠==--=nj j j j nj nj i i i j in y x l y x x x x x L 000)()(中的)(x L n 改写成:))...((...))(()()(10102010---++--=-+=n n n x x x x a x x x x a x x a a x N ,其中,n a a a ⋅⋅⋅,,10为待定定系数。
又)(00x f a =,[]0001)()(,x x x x f x f x x f a ---==,[][][]nn n n n x x x x x f x x x f x x x f a --==-,...,,,...,,,...,,10100。
将n a a a ⋅⋅⋅,,10带入)(x N n 可得:))...()(](,...,,[...)](,[)()(110100100----++-+=n n x x x x x x x x x f x x x x f x f x f 。
三、Matlab 程序代码:function newton_interpolation(x0,x) format long n=length(x0); syms zf=sqrt(1+cosh(z)^2); a(1)=subs(f,z,x0(1)); for k=1:n-1y0=subs(f,z,x0(k)); y1=subs(f,z,x0(k+1));d(k,1)=(y1-y0)/(x0(k+1)-x0(k));%一阶差商 endfor j=2:n-1for k=1:n-jd(k,j)=(d(k+1,j-1)-d(k,j-1))/(x0(k+j)-x0(k));%二阶差商及以上endenddouble(d)for j=2:na(j)=d(1,j-1);endb(1)=1;c(1)=a(1);for j=2:nb(j)=(x-x0(j-1)).*b(j-1);c(j)=a(j).*b(j);endnp=double(sum(c))wucha=double(abs(np-subs(f,z,x)))四、运行结果:五、Newton插值图像绘制实验六:一、实验名称及题目:Gauss 求积公式例5.7(P140):试构造Gauss 型求积公式∫f(x)dx 1−1≈A 0f (x 0)+A 1f (x 1)+A 2f (x 2), 并由此计算积分∫√t (1+t)2dt 10. 二、解题思路:设高斯-勒让德求积公式是:()⎪⎪⎭⎫ ⎝⎛++⎪⎪⎭⎫ ⎝⎛-≈⎰53053)(21011-f A f A f A dx x f ,依次代入2,,1)(x x x f =,解得88,95120===A A A 。