当前位置：文档之家 › c++数值分析上机实验

c++数值分析上机实验

  • 格式：doc
  • 大小：196.00 KB
  • 文档页数：24

张体强 1026222 海测一第一章舍入误差和数值稳定性1：例题利用递推公式计算定积分⎰+1n5xdxx方法一：利用递归公式151--=n n y ny （n=0,1,2,3…….20）方法二：利用递推公式n n y ny 51511-=-（n=0,1,2,3…….20）这两种方法的思路是： 1：逆代初值2：根据范围来确定循环次数 3：根据逆代函数确定逆代形似 4：如何排版1.1：方法一的C 语言实现： #include <stdio.h> #include<conio.h> #include<math.h> main() {float y0=log(6.0)-log(5.0),y1;int n=1;printf("y[0]=%-20f",y0);while(n<=20){y1=1.0/n-5*y0;printf("y[%d]=%-20f",n,y1);y0=y1;n++;if(n%3==0)printf("\n");}}1.2：方法二的c语言实现：#include <stdio.h>#include<math.h>main(){float y20=(1/105.0+1/126.0)/20,y19;int n=20;printf("y[20]=%-20f",y20);while(n>0){y19=1.0/(5*n)-y20/5;printf("y[%d]=%-20f",n-1,y19);y20=y19;n--;if(n%3==0)printf("\n");}}通过两种计算方式可以看出第二种的算法收敛2：课后习题13根据pn=(10/3)*p n-1-p n-2可知道pn的相对误差限为en=(10/3)*e n-1-e n-2 n=2 3 …….根据书中的例题方法很难比较误差是否增大，现在经过以下递归程序可以看出其结果：#include <stdio.h>double fact(int n);void main(){int n;double result;printf("Input n:");scanf("%d",&n);result=fact(n);if (result==-1)printf("Input error");/* n的值输入不正确*/elseprintf("wc=%.5f\n",result);}double fact(int n){if(n<2)return -1;elseif(n==2) /* 终止条件*/return 0.005;/* 起始误差为0.005*/elsereturn (10.0/3)*fact(n-1)-fact(n-2);/* 调用递归*/}结果表明当n=4时误差要大于n=1时的误差，所以不稳定。

17：编写程序如下#include <stdio.h>void main(){float i,j,n;float sum1=0;float sum2=0;float zhi,wc1,wc2;for(i=2;i<=100;i++){sum1+=1/(i*i-1);}printf("sum1=%f\n",sum1);for(j=100;j>=2;j--){sum2+=1/(j*j-1);}printf("sum2=%f\n",sum2);zhi=1.0/2*(3.0/2-(1.0/n)-1.0/(n+1));printf("zhi=%f\n",zhi);wc1=zhi-sum1;wc2=zhi-sum2;printf("wc1=%f\n",wc1);printf("wc2=%f\n",wc2);}三种情况结果对比：通过上图可以看出：按从大到小的顺序计算的值与精确值有较大的误差，而按从小到大的顺序计算的值与精确值较吻合，误差较小。

计算机在进行数值计算时会出现“大数吃小数”的现象，导致计算结果的精度有所降低，我们在计算机中进行同号数的加法时，采用绝对值较小者先加的算法，其结果的相对误差较小。

第二章二分法查找和牛顿逆代1：二分法查找算法设计给定区间[a b]，并设f(a)与f(b)符号相反，取ε为根的容许误差，δ为|f(x)|的容许误差1：令c=(a+b)/22：如果（c - a）<ε或就输出结果，否者执行33：如果f(a)f(c)>0,则把a的值和函数值复制给c和c的函数值，否者就把b的值和函数值复制给c和c的函数值1.1：例题：求方程分f(x)=X^3+X^2-3X-3=0在1.5附近的根（以后在这上面修改一下即可）#include <stdio.h>#include<math.h>#include<conio.h>#include<stdlib.h>#include <windows.h>#define EPS 5e-6#define wucha 1e-6float erfenfa(float,float,float(*f)(float));float f(float x);void main(){float a=1,b=2;float x;system("cls");//清除文本模式窗口x=erfenfa(a,b,f);printf("\n原方程的根是：%f",x);getch();//：从控制台读取一个字符，但不显示在屏幕上}float erfenfa(float a,float b,float(*f)(float)){float c,fc,fa=f(a),fb=f(b);int n=1;printf("二分法查找次数c\t\t f(c)\n");while(1){if(fa*fb>0) exit(1);c=(a+b)/2,fc=f(c);if(fabs(fc)<wucha) break;else if(fa*fc<0){b=c;fb=fc;}else{a=c;fa=fc;}if(b-a<EPS)break;printf("%d \t\t%f\t\t%f\n",n++,c,fc);}return c;}float f(float x){return x*x*x+x*x-3*x-3;}2：牛顿逆代算法设计给定初始值x0,ε为根的容许误差，η为|f（x）|的容许误差，n为逆代次数的容许值。

1 计算x1=x0-f(x0)/f’(x0)。

2 if(eps(x1-x0)<ε则输出x1，程序结束，否者执行33 令x0=x1,回头执行1利用c语言和matlab实现第二十题20：（1）：方程为x*x*x/3-1,误差靠输入//在纸上画出图形，估计三个值的范围//编写一个程序，用牛顿逆代求根.#include "stdafx.h"#include<stdio.h>#include<math.h>float fx(float x);float fdx(float x);void main(){float x0,x1,wc,m;int k=0;printf("请输入允许误差");/*输入逆代误差*/scanf("%f",&wc);printf("请输入初值,初值不能为-1，1\n");/*初值为1和-1导数没有意义*/ scanf("%f",&x0);do{m=-fx(x0)/fdx(x0);x1=x0+m;x0=x1;}while(fabs(m)>wc);printf("逆代结果为:%f",x0);}float fx(float x)/*函数表达式*/{float y;y=x*x*x/3-x;return(y);}float fdx(float x)/*导数*/{float yd;yd=x*x-1;return(yd);}运行结果：一个根为：1.73205一个根为：-1.73205一个根为：0（2）题目的意思和方法待请教老师。

利用matlab编写（1）如下x=-4:0.1:4;y=x.^3./3-x;plot(x,y,'r');grid on;function x=newton(fname,dfname,x0,e)if nargin<4,e=1e-4;endx=x0;x0=x+2*e;while abs(x0-x)>ex0=x;x=x0-feval(fname,x0)/feval(dfname,x0);end%定义牛顿逆代fun=inline('x^3/3-x');dfun=inline('x^2-1');format long;a=newton(fun,dfun,-2,1e-6),format short %当在-2附近时，求根a =-1.73205080756888%复制上面代码，修改求其当在0附近时，求根fun=inline('x^3/3-x');dfun=inline('x^2-1');format long;b=newton(fun,dfun,0,1e-6),format shortb=%复制上面代码，修改求其在2附近时，求根fun=inline('x^3/3-x');dfun=inline('x^2-1');format long;c=newton(fun,dfun,2,1e-6),format short c =1.73205080756888利用matlab ，在定义好牛顿逆代后，只需要三句代码可求其根，方法简单，最主要的是牛顿逆代可以以后再次调用第三章线性方程组的数值解法（利用c 语言和matlab 实现）3.1 一般方法及AX=B 得到X=A -1B例如方程 104325848332=++=++=++z y x zy x z y x 的解 输入代码如下A=[1 2 3;8 4 8;1 3 4]; B=[3;25;10]; X=A\B 结果为 X = 10.2500 28.2500 -21.2500所以x=10.2500, y=28.2500, z=-21.25003.2 列主元高斯消去法（matlab 和c 语言的实现） 3.2.1 Matlab 实现列主元高斯消去104325848332=++=++=++z y x z y x z y x在matlab 输入一下程序即可（实现编号gauss 保存在自己的常用函数盘里，这里代码省去）A=[1 2 3;8 4 8;1 3 4]; B=[3;25;10]; gauss(A,B) 结果为 x=10.2500 y=28.2500 z=-21.25003.2.2 c 语言实现列主元高斯消去求解（可参考《计算方法》3版）例如方程 104325848332=++=++=++z y x z y x z y x 的解 输入代码如下1 将方程的增广矩阵[A|b]=(a ij )n x (n+1)表示2 消元过程 对 k=1,2,3,4…………n-13 选主元4 交换5 消元#include<stdlib.h>#include <stdio.h>#include<conio.h>#include<malloc.h>#include<math.h>void main(){int i;float *x;float c[3][4]={1,2,3,3,8,4,8,25,1,3,4,10};float *liezhuyuan(float *c,int n);x=liezhuyuan(c[0],3);for(i=0;i<=2;i++)printf("x[%d]=%f\n",i,x[i]);getch();}float *liezhuyuan(float *c,int n){int i,j,t,k;float *x,p;x=(float *)malloc(n*sizeof(float));for(i=0;i<=n-2;i++){k=i;for(j=j+1;j<=n-1;j++);if(fabs(*(c+j*(n+1)+i))>(fabs(*(c+k*(n+1)+i)))) k=j;if(k==i)for(j=i;j<=n;j++){p=*(c+i*(n+1)+j);*(c+i*(n+1)+j)=*(c+k*(n+1)+j);*(c+k*(n+1)+j)=p;}for(j=i+1;j<=n-1;j++){p=(*(c+j*(n+1)+i))/(*(c+i*(n+1)+i));for(t=i;t<=n;t++)*(c+j*(n+1)+t)-=p*(*(c+i*(n+1)+t));}}for(i=n-1;i>=0;i--){for(j=n-1;j>=i+1;j--)(*(c+i*(n+1)+n))-=x[j]*(*(c+i*(n+1)+j)); x[i]=*(c+i*(n+1)+n)/(*(c+i*(n+1)+i)); } return x; }3.3 雅可比逆代（matlab 代码略）例如方程 104325848332=++=++=++z y x z y x z y x 的解 输入代码如下#include<stdio.h> #include<conio.h>#include<malloc.h> #include<math.h> #define EPS 1e-6 #define MAX 100float *jacobi(float a[3][4],int n)//定义了一个指针函数 {float *x,*y,epshe,s; int i,j,k=0;x=(float *)malloc(n*sizeof(float));//说明有几个未知数y=(float *)malloc(n*sizeof(float));//逆代后结果for(i=0;i<n;i++)x[i]=0;while(1){epshe=0;k++;for(i=0;i<n;i++){s=0;for(j=0;j<n;j++){if(j==i)continue;s+=a[i][j]*x[j];}y[i]=(a[i][n]-s)/a[i][i];epshe+=fabs(y[i]-x[i]);//y[i]为逆代值，x[i]为前面的逆代得结果}if(epshe<EPS){printf("逆代次数为：%d\n",k);return y;}if(k>=MAX){printf("这个逆代结果为发散");return y;}for(i=0;i<n;i++)x[i]=y[i];//把逆代的值赋给y[i]，以便后面的逆代，直到输出}}main(){int i;float a[3][4]={1,2,3,3,8,4,8,25,1,3,4,10};float *x;x=(float *)malloc(3*sizeof(float));x=jacobi(a,3);for(i=0;i<3;i++)printf("\nx[%d]=%f\n",i,x[i]);getch();}所以上面的方程不能雅可比逆代格式，当然有的方程可以例如：修改一下方程组的系数float a[3][4]={5,2,1,8,2,8,-3,21,1,-3,-6,1};结果为第四章插值与最小二乘发拟合牛顿逆代和Lagrange结果一样，这里我只研究Lagrange逆代。

合集下载

数值分析2024上机实验报告

数值分析2024上机实验报告

数值分析2024上机实验报告数值分析是计算数学的一个重要分支，它研究如何用数值方法来解决数学问题。

在数值分析的学习过程中，学生需要通过上机实验来巩固理论知识，并学会使用相应的数值方法来解决实际问题。

本篇报告将详细介绍2024年度数值分析上机实验的内容和结果。

一、实验内容2024年度数值分析上机实验分为四个部分，分别是：方程求根、插值与拟合、数值积分和常微分方程的数值解。

1.方程求根这部分实验要求使用数值方法求解给定的非线性方程的根。

常见的数值方法有二分法、牛顿法、割线法等。

在实验过程中，我们需要熟悉这些数值方法的原理和实现步骤，并对不同方法的收敛性进行分析和比较。

2.插值与拟合这部分实验要求使用插值和拟合方法对给定的一组数据进行拟合。

插值方法包括拉格朗日插值、牛顿插值等；拟合方法包括最小二乘拟合、多项式拟合等。

在实验中，我们需要熟悉插值和拟合方法的原理和实现步骤，并对不同方法的精度和稳定性进行比较。

3.数值积分这部分实验要求使用数值方法计算给定函数的积分。

常见的数值积分方法有梯形法则、辛普森法则、龙贝格积分等。

在实验过程中，我们需要熟悉这些数值积分方法的原理和实现步骤，并对不同方法的精度和效率进行比较。

4.常微分方程的数值解这部分实验要求使用数值方法求解给定的常微分方程初值问题。

常见的数值方法有欧拉法、改进的欧拉法、四阶龙格-库塔法等。

在实验中，我们需要熟悉这些数值解方法的原理和实现步骤，并对不同方法的精度和稳定性进行比较。

二、实验结果在完成2024年度数值分析上机实验后，我们得到了以下实验结果：1.方程求根我们实现了二分法、牛顿法和割线法，并对比了它们的收敛速度和稳定性。

结果表明，割线法的收敛速度最快，但在一些情况下可能会出现振荡；二分法和牛顿法的收敛速度相对较慢，但稳定性较好。

2.插值与拟合我们实现了拉格朗日插值和最小二乘拟合，并对比了它们的拟合效果和精度。

结果表明，拉格朗日插值在小区间上拟合效果较好，但在大区间上可能出现振荡；最小二乘拟合在整体上拟合效果较好，但可能出现过拟合。

数值分析上机实习报告

数值分析上机实习报告

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

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

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

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

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

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

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

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

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

东南大学数值分析上机题C语言编程源程序(前三章)

东南大学数值分析上机题C语言编程源程序(前三章)

第一章作业#include<stdio.h>const int N=100;int main(){int n;float S1,S2;double S3;S1=0;S2=0;for(n=2;n<N+1;n++){S1=S1+1.0/(n*n-1);}printf("the sum S1 is %d\n",S1);for(n=N;n>1;n--){S2=S2+1.0/(n*n-1);}printf("the sum S2 is %d\n",S2);S3=0.5*(1.5-1.0/N-1.0/(N+1));printf("S3 is %d\n",S3);}修改N的值即可第二章作业（1）float df(float x){float df;df=x*x-1;return df;}main(){float x0,x1,a;int k=0;printf("请输入初值x0=");scanf_s("%f",&x0);do{a=-f(x0)/df(x0);x1=x0+a;k++;x0=x1;}while(fabs(a)>eps);printf("迭代次数k=%d,根值x0=%f",k,x0);}（2）#include<stdio.h>#include<math.h>#define eps 0.01float f(float x){float f;f=x*x*x/3-x;return f;}float df(float x){float df;df=x*x-1;return df;}float ddf(float x){float ddf;ddf=2*x;return ddf;}main(){float x0,x1,a,dt=0;int k=0;do{dt=dt+eps;k++;}while((f(-dt)*f(dt)<0)&&(df(dt)!=0)&&(2*dt>=-f(-dt)/df(-dt))&&(2*dt>=f(dt)/df(dt)));printf("请输入x0：");scanf_s("%f",&x0);if(x0<-1)do{a=-f(x0)/df(x0);x1=x0+a;k++;x0=x1;}while(fabs(a)>eps);else if(x0>-1&&x0<-dt)do{a=-f(x0)/df(x0);x1=x0+a;k++;x0=x1;}while(fabs(a)>eps);else if(x0>-dt&&x0<dt)do{a=-f(x0)/df(x0);x1=x0+a;k++;x0=x1;}while(fabs(a)>eps);else if(x0>dt&&x0<1)do{a=-f(x0)/df(x0);x1=x0+a;k++;x0=x1;}while(fabs(a)>eps);else if (x0>1)do{a=-f(x0)/df(x0);x1=x0+a;k++;x0=x1;}while(fabs(a)>eps);printf("收敛域=%f,迭代次数=%d\n",dt,k);printf("根=%f\n",x0);}第三章作业#include<stdio.h>#include<math.h>main(){int n,i,j,k,t;float a[10][11],s[10],x[10],s2,max,sum,b;printf("请输入n= ");scanf("%d",&n);printf("请输入矩阵A= \n");for(i=0;i<=n;i++)for(j=0;j<=(n+1);j++)scanf("%f",&a[i][j]);//变换成高斯消去法矩阵//for(k=0;k<n;k++){//找出第k列中绝对值最大的一项//j=0;for(i=k;i<=n;i++){s[j]=a[i][k];j++;}max=fabs(s[0]);t=k;for(i=1;i<j;i++){if (fabs(s[i])>max){max=fabs(s[i]);t=k+i;}}// printf("t=%d\n",t);//将在列中有最大的元素的一行和第k行交换//for(j=0;j<=n+1;j++){b=a[k][j];a[k][j]=a[t][j];a[t][j]=b;}/*printf("***\n");for(i=0; i<=n; i++){for(j=0; j<=(n+1); j++)printf("%f ",a[i][j]);printf("\n");}*///将第k行k列的元素变为零//for(i=(k+1);i<=n;i++){s2=a[i][k]/a[k][k];for(j=k;j<=(n+1);j++)a[i][j]=a[i][j]-s2*a[k][j];}/*printf("###\n");for(i=0; i<=n; i++){for(j=0; j<=(n+1); j++)printf("%f ",a[i][j]);printf("\n");}*/printf("\n");}printf("\n");//解方程--//for(i=n;i>=0;i--){sum=0;x[n+1]=a[n][n+1]/a[n][n];for(j=i;j<=n;j++)sum+=a[i-1][j]*x[j+1];x[i]=(a[i-1][n+1]-sum)/a[i-1][i-1];}printf("方程的解为：");for(i=0;i<=n;i++)printf("%f,",x[i+1]);}。

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

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

数值分析上机实验报告《数值分析》上机实验报告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.使用此方法求方解，用误差来控制循环迭代次数，可以在误差允许的范围内得到比较理想的计算结果。