东南大学《数值分析》上机题

合集下载

东南大学数值分析上机题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]);}。

数值分析上机题目

数值分析上机题目4(总21页) --本页仅作为文档封面，使用时请直接删除即可----内页可以根据需求调整合适字体及大小--实验一实验项目：共轭梯度法求解对称正定的线性方程组实验内容：用共轭梯度法求解下面方程组(1) 123421003131020141100155x x x x -⎛⎫⎛⎫⎛⎫ ⎪ ⎪ ⎪--- ⎪ ⎪ ⎪=⎪ ⎪ ⎪-- ⎪ ⎪ ⎪-⎝⎭⎝⎭⎝⎭ 迭代20次或满足()(1)1110k k x x --∞-<时停止计算。

编制程序：储存m 文件function [x,k]=CGmethod(A,b)n=length(A);x=2*ones(n,1);r=b-A*x;rho=r'*r; k=0;while rho>10^(-11) & k<1000 k=k+1; if k==1 p=r; elsebeta=rho/rho1; p=r+beta*p; end w=A*p;alpha=rho/(p'*w); x=x+alpha*p; r=r-alpha*w; rho1=rho;rho=r'*r; end运行程序： clear,clcA=[2 -1 0 0;-1 3 -1 0;0 -1 4 -1;0 0 -1 5]; b=[3 -2 1 5]'; [x,k]=CGmethod(A,b)运行结果： 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 -=-≤时停止计算。

编制程序：储存m 文件function [x,k]=CGmethod_1(A,b) n=length(A);x(1:n,1)=0;r=b-A*x;r1=r; k=0;while norm(r1,1)>10^(-7)&k<10^4 k=k+1; if k==1 p=r; elsebeta=(r1'*r1)/(r'*r);p=r1+beta*p; end r=r1; w=A*p;alpha=(r'*r)/(p'*w); x=x+alpha*p; r1=r-alpha*w; end运行程序： clear,clc n=1000; A=hilb(n); b=sum(A')';[x,k]=CGmethod(A,b)实验二1、实验目的：用复化Simpson 方法、自适应复化梯形方法和Romberg 方法求数值积分。

数值分析报告上机题课后作业全部-东南大学

实用标准文案文档大全上机作业题报告2015.1.9 USER1.Chapter 11.1题目设S N =∑1j 2−1N j=2，其精确值为)11123(21+--N N 。

（1）编制按从大到小的顺序11131121222-+⋯⋯+-+-=N S N ，计算S N 的通用程序。

（2）编制按从小到大的顺序1211)1(111222-+⋯⋯+--+-=N N S N ，计算S N 的通用程序。

（3）按两种顺序分别计算64210,10,10S S S ,并指出有效位数。

（编制程序时用单精度）（4）通过本次上机题，你明白了什么？1.2程序1.3运行结果1.4结果分析按从大到小的顺序，有效位数分别为：6，4，3。

按从小到大的顺序，有效位数分别为：5，6，6。

可以看出，不同的算法造成的误差限是不同的，好的算法可以让结果更加精确。

当采用从大到小的顺序累加的算法时，误差限随着N 的增大而增大，可见在累加的过程中，误差在放大，造成结果的误差较大。

因此，采取从小到大的顺序累加得到的结果更加精确。

2.Chapter 22.1题目（1）给定初值0x 及容许误差ε，编制牛顿法解方程f(x)=0的通用程序。

（2）给定方程03)(3=-=x xx f ,易知其有三个根3,0,3321=*=*-=*x x x○1由牛顿方法的局部收敛性可知存在,0>δ当),(0δδ+-∈x 时，Newton 迭代序列收敛于根x2*。

试确定尽可能大的δ。

○2试取若干初始值，观察当),1(),1,(),,(),,1(),1,(0+∞+-----∞∈δδδδx 时Newton 序列的收敛性以及收敛于哪一个根。

（3）通过本上机题，你明白了什么？2.2程序2.3运行结果（1）寻找最大的δ值。

算法为：将初值x0在从0开始不断累加搜索精度eps，带入Newton迭代公式，直到求得的根不再收敛于0为止，此时的x0值即为最大的sigma值。

运行Find.m，得到在不同的搜索精度下的最大sigma值。

东南大学数值分析上机实验题(下)

数值分析上机报告XX：学号：2013年12月22日第四章38.（上机题）3次样条插值函数（1）编制求第一型3次样条插值函数的通用程序;端点条件为'0y =0.8，'10y =0.2。

用所编制程序求车门的3次样条插值函数S(x),并打印出S(i+0.5)(i=0,1,…9)。

解：（1）#include <iostream.h> #include <math.h>floatx1[100],f1[100],f2[99],f3[98],m[100],a[100][101],x,d[100]; float c[99],e[99],h[99],u[99],w[99],y_0,y_n,arr ,s; int i,j,k,n,q;void selectprint(float y) {if ((y>0)&&(y!=1)) cout<<"+"<<y; else if (y==1) cout<<"+"; else if (y<0) cout<<y; }void printY(float y){ if (y!=0) cout<<y; }float calculation(float x){ for (j=1;j<=n;j++) if (x<=x1[j]) {s=(float)(f1[j-1]+c[j-1]*(x-x1[j-1])+m[j-1]/2.0*(x-x1[j-1])*(x-x1[j-1])+e[j-1]*(x-x1[j-1])*(x-x1[j-1])*(x-x1[j-1])); break; }return s; }void main() {do{cout<<"请输入n值:";cin>>n;if ((n>100)||(n<1)) cout<<"请重新输入整数(1..100):"<<endl;}while ((n>100)||(n<1));cout<<"请输入Xi(i=0,1,...,"<<n<<"):";for (i=0;i<=n;i++) cin>>x1[i];cout<<endl;cout<<"请输入Yi(i=0,1,...,"<<n<<"n):";for (i=0;i<=n;i++) cin>>f1[i];cout<<endl;cout<<"请输入第一型边界条件Y'0，Y'n:";cin>>y_0>>y_n;cout<<endl;for (i=0;i<n;i++) h[i]=x1[i+1]-x1[i];for (i=1;i<n;i++) u[i]=(float) (h[i-1]/(h[i-1]+h[i]));for (i=1;i<n;i++) w[i]=(float) (1.0-u[i]);for (i=0;i<n;i++) f2[i]=(f1[i+1]-f1[i])/h[i]; //一阶差商for (i=0;i<n-1;i++) f3[i]=(f2[i+1]-f2[i])/(x1[i+2]-x1[i]); //二阶差商for (i=1;i<n;i++) d[i]=6*f3[i-1]; //求出所有的d[i]d[0]=6*(f2[0]-y_0)/h[0];d[n]=6*(y_n-f2[n-1])/h[n-1];for (i=0;i<=n;i++)for (j=0;j<=n;j++)if (i==j) a[i][j]=2;else a[i][j]=0;a[0][1]=1;a[n][n-1]=1;for (i=1;i<n;i++){a[i][i-1]=u[i];a[i][i+1]=w[i];}for (i=0;i<=n;i++) a[i][n+1]=d[i];for (i=1;i<=n;i++) //用追赶法解方程,得M[i]{arr=a[i][i-1];for (j=0;j<=n+1;j++)a[i][j]=a[i][j]-a[i-1][j]*arr/a[i-1][i-1];}m[n]=a[n][n+1]/a[n][n];for (i=n-1;i>=0;i--) m[i]=(a[i][n+1]-a[i][i+1]*m[i+1])/a[i][i];for (i=0;i<n;i++) //计算S（x）的表达式c[i]=(float) (f2[i]-(1.0/3.0*m[i]+1.0/6.0*m[i+1])*h[i]);for (i=0;i<n;i++)e[i]=(m[i+1]-m[i])/(6*h[i]);for (i=0;i<n;i++){cout<<"X属于区间["<<x1[i]<<","<<x1[i+1]<<"]时"<<endl<<endl;cout<<"S(x)=";printY(f1[i]);if (c[i]!=0){selectprint(c[i]);cout<<"(x";printY(-x1[i]);cout<<")";}if (m[i]!=0){selectprint((float)(m[i]/2.0));for (q=0;q<2;q++){cout<<"(x";printY(-x1[i]);cout<<")";}}if (e[i]!=0){selectprint(e[i]);for (q=0;q<3;q++){cout<<"(x";printY(-x1[i]);cout<<")";}}cout<<endl<<endl;}cout<<"针对本题，计算S(i+0.5)(i=0,1,...,9):"<<endl;for (i=0;i<10;i++){if ((i+0.5<=x1[n])&&(i+0.5>=x1[0])){calculation((float)(i+0.5));cout<<"S("<<i+0.5<<")="<<s<<endl;}else cout<<i+0.5<<"超出定义域"<<endl;};cout<<endl;}（2）编制的程序求车门的3次样条插值函数S(x)：x属于区间[0,1]时；S(x)=2.51+0.8(x)-0.0014861(x)(x)-0.00851395(x)(x)(x)x属于区间[1,2]时；S(x)=3.3+0.771486(x-1)-0.027028(x-1)(x-1)-0.00445799(x-1)(x-1)(x-1) x属于区间[2,3]时；S(x)=4.04+0.704056(x-2)-0.0404019(x-2)(x-2)-0.0036543(x-2)(x-2)(x-2)x属于区间[3,4]时；S(x)=4.7+0.612289(x-3)-0.0513648(x-3)(x-3)-0.0409245(x-3)(x-3)(x-3) x属于区间[4,5]时；S(x)=5.22+0.386786(x-4)-0.174138(x-4)(x-4)+0.107352(x-4)(x-4)(x-4) x属于区间[5,6]时；S(x)=5.54+0.360567(x-5)+0.147919(x-5)(x-5)-0.268485(x-5)(x-5)(x-5) x属于区间[6,7]时；S(x)=5.78-0.149051(x-6)-0.657537(x-6)(x-6)+0.426588(x-6)(x-6)(x-6) x属于区间[7,8]时；S(x)=5.4-0.184361(x-7)+0.622227(x-7)(x-7)-0.267865(x-7)(x-7)(x-7) x属于区间[8,9]时；S(x)=5.57+0.256496(x-8)-0.181369(x-8)(x-8)+0.0548728(x-8)(x-8)(x-8) x属于区间[9,10]时；S(x)=5.7+0.058376(x-9)-0.0167508(x-9)(x-9)+0.0583752(x-9)(x-9)(x-9)S(0.5)=2.90856 S(1.5)=3.67843 S (2.5)=4.38147 S(3.5)=4.98819 S(4.5)=5.38328 S(5.5)=5.7237S(6.5)=5.59441 S(7.5)=5.42989 S(8.5)=5.65976 S(9.5)=5.7323第六章21．（上机题）常微分方程初值问题数值解（1）编制RK4方法的通用程序；（2）编制AB4方法的通用程序（由RK4提供初值）；（3）编制AB4-AM4预测校正方法的通用程序（由RK4提供初值）；（4）编制带改进的AB4-AM4预测校正方法的通用程序（由RK4提供初值）;（5）对于初值问题h=,应用（1）~（4）中的四种方法进行计算，并将计算结果和精确解取步长0.13=+作比较；y x x()3/(1)（6）通过本上机题，你能得到哪些结论？解：#include<iostream.h>#include<fstream.h>#include<stdlib.h>#include<math.h>ofstream outfile("data.txt");//此处定义函数f(x,y)的表达式//用户可以自己设定所需要求得函数表达式double f1(double x,double y){double f1;f1=(-1)*x*x*y*y;return f1;}//此处定义求函数精确解的函数表达式double f2(double x){double f2;f2=3/(1+x*x*x);return f2;}//此处为精确求函数解的通用程序void accurate(double a,double b,double h){double x[100],accurate[100];x[0]=a;int i=0;outfile<<"输出函数准确值的程序结果:\n";do{x[i]=x[0]+i*h;accurate[i]=f2(x[i]);outfile<<"accurate["<<i<<"]="<<accurate[i]<<'\n';i++;}while(i<(b-a)/h+1);}//此处为经典Runge-Kutta公式的通用程序void RK4(double a,double b,double h,double c) {int i=0;double k1,k2,k3,k4;double x[100],y[100];y[0]=c;x[0]=a;outfile<<"输出经典Runge-Kutta公式的程序结果:\n"; do{x[i]=x[0]+i*h;k1=f1(x[i],y[i]);k2=f1((x[i]+h/2),(y[i]+h*k1/2));k3=f1((x[i]+h/2),(y[i]+h*k2/2));k4=f1((x[i]+h),(y[i]+h*k3));y[i+1]=y[i]+h*(k1+2*k2+2*k3+k4)/6;outfile<<"y"<<"["<<i<<"]="<<y[i]<<'\n';i++;}while(i<(b-a)/h+1);}//此处为4阶Adams显式方法的通用程序void AB4(double a,double b,double h,double c) {double x[100],y[100],y1[100];double k1,k2,k3,k4;y[0]=c;x[0]=a;outfile<<"输出4阶Adams显式方法的程序结果:\n";for(int i=0;i<=2;i++){x[i]=x[0]+i*h;k1=f1(x[i],y[i]);k2=f1((x[i]+h/2),(y[i]+h*k1/2));k3=f1((x[i]+h/2),(y[i]+h*k2/2));k4=f1((x[i]+h),(y[i]+h*k3));y[i+1]=y[i]+h*(k1+2*k2+2*k3+k4)/6;}int j=3;y1[0]=y[0];y1[1]=y[1];y1[2]=y[2];y1[3]=y[3];do{x[j]=x[0]+j*h;y1[j+1]=y1[j]+(55*f1(x[j],y1[j])-59*f1(x[j-1],y1[j-1])+37*f1(x[j-2],y1[j-2])-9*f1(x[j-3],y1[ j-3]))*h/24;outfile<<"y1"<<"["<<j<<"]="<<y1[j]<<'\n';j++;}while(j<(b-a)/h+1);}//主函数void main(void){double a,b,h,c;cout<<"输入上下区间、步长和初始值:\n";cin>>a>>b>>h>>c;accurate(a,b,h);RK4(a,b,h,c);AB4(a,b,h,c);}float si(int u,int v){float sum=0; int q;for(q=0;q<k;q++)sum+=a[u][q]*a[q][v];sum=a[u][v]-sum;return sum;}void exchange(int g){int t; float temp;for(t=0;t<n;t++){temp=a[k][t];a[k][t]=a[g][t];a[g][t]=temp;}}void analyze(){int t;float si(int u,int v);for(t=k;t<n;t++)a[k][t]=si(k,t);for(t=(k+1);t<m;t++)a[t][k]=(float)(si(t,k)/a[k][k]);}void ret(){int t,z;float sum;x[m-1]=(float)a[m-1][m]/a[m-1][m-1];for(t=(m-2);t>-1;t--){sum=0;for(z=(t+1);z<m;z++)sum+=a[t][z]*x[z];x[t]=(float)(a[t][m]-sum)/a[t][t];}}（5）由经典Runge-Kutta公式得出的结果列在下面的表格中，以及精确值y(x i)和精确值和数值解的误差：由AB4方法得出的结果为：Y1[0]=3 y1[1]=2.997 y1[2]=2.97619 y1[3]=2.92113 y1[4]=2.81839 y1[5]=2.66467 y1[6]=2.4652 y1[7]=2.23308 y1[8]=1.98495 y1[9]=1.73704 y1[10]=1.50219 y1[11]=1.28876 y1[12]=1.10072 y1[13]=0.93871 y1[14]=0.801135y1[15]=0.685335（6）本次上机作业让我知道了在遇到复杂问题中，无法给出解析式的情况下，要学会灵活使用各种微分数值解法，并且能计算出不同方法的精度大小。

东南大学数值分析上机练习后三章

数值分析上机练习（以VC++6.0为操作平台）第四章（4.38）程序如下：#include<iostream.h>void main(void){float x[11];//存放数组x[j]float y[11];//存放数组y[j]float h[11];//存放数组h[j]float u[11];//存放数组u[j]float v[11];//存放数组v[j]float d[11];//存放数组d[j]float M[11];//存放数组M[j]float b[11];// 存放数组b[j]float t[11],l[11],yy[11],s[4],aa1,aa2,aa3,aa4;float s1[10];int i,j,n;float xx;cout<<"请输入n的值:\n";cin>>n;cout<<"输入数组x:\n";for(i=0;i<=n;i++)cin>>x[i];cout<<"输入数组y:\n";for(i=0;i<=n;i++)cin>>y[i];//输入端点值float df[2];cout<<"输入两个端点值:\n";for(i=0;i<2;i++)cin>>df[i];//求出h[j]的值for(j=0;j<=n-1;j++){h[j]=x[j+1]-x[j];cout<<'h'<<'['<<j<<']'<<'='<<h[j]<<'\t';}cout<<endl;//求出u[j]和v[j]的初值v[0]=1;u[n]=1;for(j=1;j<=n-1;j++){u[j]=h[j-1]/(h[j-1]+h[j]);v[j]=h[j]/(h[j-1]+h[j]);}//求出d[j]的值for(j=1;j<n;j++){d[j]=6*((y[j+1]-y[j])/h[j]-(y[j]-y[j-1])/h[j-1])/(h[j]+h[j-1]);} d[0]=6*((y[1]-y[0])/h[0]-df[0])/h[0];d[n]=6*(df[1]-(y[n]-y[n-1])/h[n-1])/h[n-1];for(j=1;j<=n;j++){cout<<'u'<<'['<<j<<']'<<'='<<u[j]<<'\t';}cout<<endl;for(j=0;j<n;j++){cout<<'v'<<'['<<j<<']'<<'='<<v[j]<<'\t';}cout<<endl;for(j=0;j<=n;j++){cout<<'d'<<'['<<j<<']'<<'='<<d[j]<<'\t';}cout<<endl;//利用书本上的追赶法求解方程组for(i=0;i<=n;i++){b[i]=2;}cout<<endl;t[0]=b[0];yy[0]=d[0];//消元过程for(i=1;i<=n;i++){l[i]=u[i]/t[i-1];t[i]=b[i]-l[i]*v[i-1];yy[i]=d[i]-l[i]*yy[i-1];}//回代过程M[n]=yy[n]/t[n];for(i=n-1;i>=0;i--){M[i]=(yy[i]-v[i]*M[i+1])/t[i];}//将M[j]的值输出for(i=0;i<=n;i++){cout<<'M'<<'['<<i<<']'<<'='<<M[i]<<endl;}//输出插值多项式的系数for(j=0;j<n;j++){s[0]=y[j];s[1]=(y[j+1]-y[j])/h[j]-(M[j]/3+M[j+1]/6)*h[j];s[2]=M[j]/2;s[3]=(M[j+1]-M[j])/(6*h[j]);cout<<"当x的值在区间"<<'x'<<'['<<j<<']'<<"到"<<'x'<<'['<<(j+1)<<']'<<"时,输出插值多项式的系数:\n";for(int k=0;k<4;k++){cout<<'s'<<'['<<k<<']'<<'='<<s[k]<<'\t';}cout<<endl;}}程序结果：详见附图4.38jpg编制的程序求车门的3次样条插值函数S(x)：x属于区间[0,1]时；S(x)=2.51+0.8(x)-0.0014861(x)(x)-0.00851395(x)(x)(x)x属于区间[1,2]时；S(x)=3.3+0.771486(x-1)-0.027028(x-1)(x-1)-0.00445799(x-1)(x-1)(x-1) x属于区间[2,3]时；S(x)=4.04+0.704056(x-2)-0.0404019(x-2)(x-2)-0.0036543(x-2)(x-2)(x-2) x属于区间[3,4]时；S(x)=4.7+0.612289(x-3)-0.0513648(x-3)(x-3)-0.0409245(x-3)(x-3)(x-3) x属于区间[4,5]时；S(x)=5.22+0.386786(x-4)-0.174138(x-4)(x-4)+0.107352(x-4)(x-4)(x-4) x属于区间[5,6]时；S(x)=5.54+0.360567(x-5)+0.147919(x-5)(x-5)-0.268485(x-5)(x-5)(x-5) x属于区间[6,7]时；S(x)=5.78-0.149051(x-6)-0.657537(x-6)(x-6)+0.426588(x-6)(x-6)(x-6) x属于区间[7,8]时；S(x)=5.4-0.184361(x-7)+0.622227(x-7)(x-7)-0.267865(x-7)(x-7)(x-7)x属于区间[8,9]时；S(x)=5.57+0.256496(x-8)-0.181369(x-8)(x-8)+0.0548728(x-8)(x-8)(x-8) x属于区间[9,10]时；S(x)=5.7+0.058376(x-9)-0.0167508(x-9)(x-9)+0.0583752(x-9)(x-9)(x-9) S(0.5)=2.90856 S(1.5)=3.67843 S (2.5)=4.38147S(3.5)=4.98819 S(4.5)=5.38328 S(5.5)=5.7237S(6.5)=5.59441 S(7.5)=5.42989 S(8.5)=5.65976S(9.5)=5.7323第六章（6.21）程序如下：#include<iostream.h>#include<fstream.h>#include<stdlib.h>#include<math.h>ofstream outfile("data.txt");//此处定义函数f(x,y)的表达式//用户可以自己设定所需要求得函数表达式double f1(double x,double y){double f1;f1=(-1)*x*x*y*y;return f1;}//此处定义求函数精确解的函数表达式double f2(double x){double f2;f2=3/(1+x*x*x);return f2;}//此处为精确求函数解的通用程序void accurate(double a,double b,double h){double x[100],accurate[100];x[0]=a;int i=0;outfile<<"输出函数准确值的程序结果:\n";do{x[i]=x[0]+i*h;accurate[i]=f2(x[i]);outfile<<"accurate["<<i<<"]="<<accurate[i]<<'\n';i++;}while(i<(b-a)/h+1);}//此处为经典Runge-Kutta公式的通用程序void RK4(double a,double b,double h,double c){int i=0;double k1,k2,k3,k4;double x[100],y[100];y[0]=c;x[0]=a;outfile<<"输出经典Runge-Kutta公式的程序结果:\n";do{x[i]=x[0]+i*h;k1=f1(x[i],y[i]);k2=f1((x[i]+h/2),(y[i]+h*k1/2));k3=f1((x[i]+h/2),(y[i]+h*k2/2));k4=f1((x[i]+h),(y[i]+h*k3));y[i+1]=y[i]+h*(k1+2*k2+2*k3+k4)/6;outfile<<"y"<<"["<<i<<"]="<<y[i]<<'\n';i++;}while(i<(b-a)/h+1);}void AB4(double a,double b,double h,double c){double x[100],y[100],y1[100];double k1,k2,k3,k4;y[0]=c;x[0]=a;outfile<<"输出4阶Adams显式方法的程序结果:\n";for(int i=0;i<=2;i++){x[i]=x[0]+i*h;k1=f1(x[i],y[i]);k2=f1((x[i]+h/2),(y[i]+h*k1/2));k3=f1((x[i]+h/2),(y[i]+h*k2/2));k4=f1((x[i]+h),(y[i]+h*k3));y[i+1]=y[i]+h*(k1+2*k2+2*k3+k4)/6;}int j=3;y1[0]=y[0];y1[1]=y[1];y1[2]=y[2];y1[3]=y[3];do{x[j]=x[0]+j*h;y1[j+1]=y1[j]+(55*f1(x[j],y1[j])-59*f1(x[j-1],y1[j-1])+37*f1(x[j-2], y1[j-2])-9*f1(x[j-3],y1[j-3]))*h/24;outfile<<"y1"<<"["<<j<<"]="<<y1[j]<<'\n';j++;}while(j<(b-a)/h+1);}//主函数void main(void){double a,b,h,c;cout<<"输入上下区间、步长和初始值:\n";cin>>a>>b>>h>>c;accurate(a,b,h);RK4(a,b,h,c);AB4(a,b,h,c);}程序结果：经典Runge-Kutta公式得出的结果列在下面的表格中，以及精确由AB4方法得出的结果为：y1[0]=3 y1[1]=2.997 y1[2]=2.97619 y1[3]=2.92113y1[4]=2.81839 y1[5]=2.66467y1[6]=2.4652 y1[7]=2.23308y1[8]=1.98495y1[9]=1.73704y1[10]=1.5021 y1[11]=1.28876y1[12]=1.10072y1[13]=0.93871y1[14]=0.801135 y1[15]=0.685335通过本上机题我明白了各种求微分方程的数值方法，经典Runge-Kutta公式，AB4方法以及AB4-AM4预测校正方法求解公式的精度是不同的。

东南大学出版社第二版《数值分析》上机作业答案(前三章)

for (i=k+1;i<N;i++) // { lik=a[i][k]/a[k][k]; //实施消去过程，得到上三角系数增广矩阵 for (j=k;j<M;j++) // { a[i][j]=a[i][j]‐lik*a[k][j]; // } } } cout<<"经列主元高斯消去法得到的数组为:"<<endl; // for (b=0;b<N;b++) // { cout<<endl; //输出经过列主元消去法处理过的系数增广矩阵 for (c=0;c<M;c++) { cout<<setw(7)<<a[b][c]; // } } cout<<endl; double x[N]; // double s; int f,g; x[N‐1]=a[N‐1][M‐1]/a[N‐1][N‐1]; // for (f=N‐2;f>=0;f‐‐) // { s=0; for (g=f+1;g<N;g++) //由上三角形的系数增广矩阵求出方程组的解 { s=s+a[f][g]*x[g]; // } x[f]=(a[f][N]‐s)/a[f][f]; // } cout<<"方程组的解为:"<<endl; for (b=0;b<N;b++) //输出方程组的解 {
1
当 n=10000 时，s3=0.7499 Press any key to continue （分析 S1 的 6 位数字中，有效位数为 4 位； S2 的所有数字都是有效数字。）当 n=1000000 时，s1=‐14.2546 当 n=1000000 时，s2=‐14.2551 当 n=1000000 时，s3=0.749999 Press any key to continue （分析： S1 的 6 位数字中，没有有效数字； S2 的 6 位数字中，没有有效数字。）由运行结果可知，当精度比较低时，按从大数开始累加到小数的计算结果的精度低于按从小数累加到大数的计算结果的精度。至于当 n=1000000 时，S1 和 S2 得出了负数结果，可能是由于循环次数过多，导致数据溢出，从而得出错误结果。习题 2 20．程序如下： //给定误差限为：0.5e‐6 //经过试算得当 delta 最大取道 0.7745966 时，迭代得到的根都收敛于 0 #include <iostream.h> #include <math.h> void main () { double x,u; int count=0; u=10.0; cout<<"请输入 x 的初值"<<endl; cin>>x; for (count=0;abs(u)>5;count++) { x=x‐(x*x*x‐3*x)/(3*(x*x‐1)); u=10000000*x; if(count>5000) { cout<<"迭代结果不收敛于 0！"<<endl; break; } } cout<<"x="<<x<<endl<<endl;

数值分析上机题Matlab(东南大学)3

0.22 0.24 0.26 0.28 0.30 0.32 0.34 0.36 0.38 0.40 0.42 0.44 0.46 0.48 0.50 0.52 0.54 0.56 0.58 0.60 0.62 0.64 0.66 0.68 0.70 0.72
152 139 128 119 110 103 96 90 85 80 76 72 68 65 62 59 56 53 51 49 47 45 43 41 39 38
========================================================================================================================
======================================================================================================================================================================== 习题 3_36 ======================================================================================================================================================================== Omega n x1 x2 x3 x4 x5 x6 x7 x8 x9
-0.71279 -0.71280 -0.71280 -0.71280 -0.71280 -0.71280 -0.71280 -0.71280 -0.71280 -0.71280 -0.71280 -0.71281 -0.71281 -0.71281 -0.71281 -0.71281 -0.71281 -0.71281 -0.71281 -0.71281 -0.71281 -0.71281 -0.71281 -0.71281 -0.71281 -0.71281

东南大学数值分析上机题答案

数值分析上机题第一章17.（上机题）舍入误差与有效数设∑=-=Nj N j S 2211，其精确值为）111-23（21+-N N 。

（1）编制按从大到小的顺序1-1···1-311-21222N S N +++=，计算N S 的通用程序；（2）编制按从小到大的顺序121···1)1(111222-++--+-=N N S N ，计算NS 的通用程序；（3）按两种顺序分别计算210S ,410S ,610S ,并指出有效位数（编制程序时用单精度）；（4）通过本上机题，你明白了什么？解：程序：（1）从大到小的顺序计算1-1···1-311-21222N S N +++=：function sn1=fromlarge(n) %从大到小计算sn1format long ; sn1=single(0); for m=2:1:nsn1=sn1+1/(m^2-1); end end（2）从小到大计算121···1)1(111222-++--+-=N N S N function sn2=fromsmall(n) %从小到大计算sn2format long ; sn2=single(0); for m=n:-1:2sn2=sn2+1/(m^2-1); end end（3）总的编程程序为： function p203()clear allformat long;n=input('please enter a number as the n:') sn=1/2*(3/2-1/n-1/(n+1));%精确值为snfprintf('精确值为%f\n',sn);sn1=fromlarge(n);fprintf('从大到小计算的值为%f\n',sn1);sn2=fromsmall(n);fprintf('从小到大计算的值为%f\n',sn2);function sn1=fromlarge(n) %从大到小计算sn1 format long;sn1=single(0);for m=2:1:nsn1=sn1+1/(m^2-1);endendfunction sn2=fromsmall(n) %从小到大计算sn2 format long;sn2=single(0);for m=n:-1:2sn2=sn2+1/(m^2-1);endendend运行结果：从而可以得到N值真值顺序值有效位数2 100.740050 从大到小0.740049 5从小到大0.740050 64 100.749900 从大到小0.749852 3从小到大0.749900 66 100.749999 从大到小0.749852 3从小到大0.749999 6（4）感想：通过本上机题，我明白了，从小到大计算数值的精确位数比较高而且与真值较为接近，而从大到小计算数值的精确位数比较低。

数值分析上机题目详解

第一章一、题目设∑=-=N N j S 2j 211，其精确值为)11123(21+--N N 。

1) 编制按从大到小的顺序11131121222-+⋯⋯+-+-=N S N ，计算S N 的通用程序。

2) 编制按从小到大的顺序1211)1(111222-+⋯⋯+--+-=N N S N ，计算S N 的通用程序。

3) 按两种顺序分别计算64210,10,10S S S ,并指出有效位数。

（编制程序时用单精度）4) 通过本次上机题，你明白了什么？二、通用程序N=input('Please Input an N (N>1):');AccurateValue=single((0-1/(N+1)-1/N+3/2)/2);Sn1=single(0);for a=2:N;Sn1=Sn1+1/(a^2-1);endSn2=single(0);for a=2:N;Sn2=Sn2+1/((N-a+2)^2-1);endfprintf('The value of Sn (N=%d)\n',N);fprintf('Accurate Calculation %f\n',AccurateValue); fprintf('Caculate from large to small %f\n',Sn1);fprintf('Caculate from small to large %f\n',Sn2);disp('____________________________________________________')三、结果从结果可以看出有效位数是6位。

感想：可以得出，算法对误差的传播有一定的影响，在计算时选一种好的算法可以使结果更为精确。

从以上的结果可以看到从大到小的顺序导致大数吃小数的现象，容易产生较大的误差，求和运算从小数到大数所得到的结果才比较准确。

数值分析上机作业

数值分析上机作业2实验一：（1） ①用不动点迭代法求()013=--=x x x f 的根，至少设计两种迭代格式，一个收敛一个发散，1210-=ε。

（2） ②对迭代格式使用Aitken 加速，观察敛散性变化。

1取递推公式31)1(+=x x ，可以得到收敛时的迭代结果为：x=(2)^(1/3); t=1; while(1)if(abs(x-(x+1)^(1/3))<10^-12) break; endx=(x+1)^(1/3);t=t+1;end t xtemp=x^3-x-1 %带回来验证下 t = 16 x =1.324717957243755 temp =-4.225952920933196e-12 加速后代码如下 x=1;x=(x*((x+1)^(1/3)+1)^(1/3)-(x+1)^(2/3))/(x-2*(x+1)^(1/3)+((x+1)^(1/3)+1)^(1/3)) t=1; while(1)if(abs(x-(x*((x+1)^(1/3)+1)^(1/3)-(x+1)^(2/3))/(x-2*(x+1)^(1/3)+((x+1)^(1/3)+1)^(1/3)))<10^-10) break; endx=(x*((x+1)^(1/3)+1)^(1/3)-(x+1)^(2/3))/(x-2*(x+1)^(1/3)+((x+1)^(1/3)+1)^(1/3)); t=t+1; if (t>100000)break; %防止运算次数过多 end endfprintf('需要%d 次',t);输出需要670次>>此处取10^10-=ε，若到10^-12次方则可能需要运行更多次取13-=x x 则迭代发散。

使用aitken 加速计算结果如下 x=2; t=1; while(1)if( abs(x-((x*((x^3-1)^3-1))-(x^3-1)^2)/(x-2*(x^3-1)+(x^3-1)^3-1))<10^-10) break; end t=t+1;x=((x*((x^3-1)^3-1))-(x^3-1)^2)/(x-2*(x^3-1)+(x^3-1)^3-1); if (t>100000)break; %防止运算次数过多 end end t xt = 108x = 1.324717956244172由此可见经过aitken 加速以后，原来发散的迭代格式收敛了。

东南大学数值分析上机

第一章一、题目设∑=-=Nj N j S 2211，其精确值为)11123(21+--N N 。

（1）编制按从大到小的顺序11131121222-+⋯⋯+-+-=N S N ，计算SN 的通用程序。

（2）编制按从小到大的顺序1211)1(111222-+⋯⋯+--+-=N N S N ，计算SN 的通用程序。

（3）按两种顺序分别计算64210,10,10S S S ,并指出有效位数。

（编制程序时用单精度）（4）通过本次上机题，你明白了什么？二、MATLAB 程序N=input('请输入N(N>1)：');AccurateValue=single((0-1/(N+1)-1/N+3/2)/2); %single 使其为单精度 Sn1=single(0); %从小到大的顺序 for a=2:N; Sn1=Sn1+1/(a^2-1); endSn2=single(0); %从大到小的顺序 for a=2:N; Sn2=Sn2+1/((N-a+2)^2-1); endfprintf('Sn 的值 (N=%d)\n',N);disp('____________________________________________________') fprintf('精确值 %f\n',AccurateValue); fprintf('从大到小计算的结果 %f\n',Sn1); fprintf('从小到大计算的结果 %f\n',Sn2);disp('____________________________________________________')三、结果请输入N(N>1)：100Sn的值(N=100)____________________________________________________精确值0.740049从大到小计算的结果0.740049从小到大计算的结果0.740050____________________________________________________请输入N(N>1)：10000Sn的值(N=10000)____________________________________________________精确值0.749900从大到小计算的结果0.749852从小到大计算的结果0.749900____________________________________________________请输入N(N>1)：1000000Sn的值(N=1000000)____________________________________________________精确值0.749999从大到小计算的结果0.749852从小到大计算的结果0.749999____________________________________________________四、结果分析可以得出，算法对误差的传播又一定的影响，在计算时选一种好的算法可以使结果更为精确。

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

一、问答题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 *︱≠21×10m-n+12、数值计算应该避免采用不稳定的算法，防止有效数字的损失. 因此，在进行数值运算算法设计过程中主要注意什么？（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*的范数很小，但该近似解仍不符合问题要求。

6、Lagrange 插值的前提条件是什么？并写出二次Lagrange 插值的基函数。

前提条件是：⎩⎨⎧≠==i j i j x j，，（01)l i .2,1,0,n j i ， = 二次Lagrange 插值的基函数：()))(())((2010210x x x x x x x x x l ----=()))(())((2101201x x x x x x x x x l ----= ()))(())((1202102x x x x x x x x x l ----=7、什么是数值积分的代数精度？如果某一个求积公式对于次数不超过m 的多项式均能准确地成立，但对于m+1次多项式就不准确成立，则称该求积公式具有m 次代数精度（或代数精确度）。

数值分析上机题答案

数值分析上机题姓名：武均学号：142648习题117．（上机题）舍入误差与有效数设2211NN j S j ==-∑，其精确值为1311221N N ⎛⎫-- ⎪+⎝⎭。

（1）编制按从大到小的顺序22211121311N S N =+++---，计算N S 的通用程序。

（2）编制按从小到大的顺序2221111(1)121N S N N =+++----，计算NS 的通用程序。

（3）按两种顺序分别计算210S ，410S ，610S ，并指出有效位数。

（编制程序时用单精度）（4）通过本上机题你明白了什么？按从大到小的顺序计算N S 的通用程序为： #include<iostream.h>float sum(float N) { float j,s,sum=0; for(j=2;j<=N;j++) { s=1/(j*j-1); sum+=s; } return sum; } 按从小到大的顺序计算N S 的通用程序为： #include<iostream.h> float sum(float N) {float j,s,sum=0; for(j=N;j>=2;j--) {s=1/(j*j-1); sum+=s; }return sum; }从大到小的顺序的值从小到大的顺序的值精确值有效位数从大到小从小到大210S 0.740049 0.74005 0.740049 6 5 410S0.749852 0.7499 0.7499 4 4 610S0.7498520.7499990.74999936通过本上机题，看出按两种不同的顺序计算的结果是不相同的，按从大到小的顺序计算的值与精确值有较大的误差，而按从小到大的顺序计算的值与精确值吻合。

从大到小的顺序计算得到的结果的有效位数少。

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

数值分析上机题答案

20．（上机题）舍入误差与有效数设∑=-=N j N j S 2211，其精确值为)11123(21+--N N 。

（1）编制按从小到大的顺序11131121222-+⋅⋅⋅+-+-=N S N ，计算N S 的通用程序。

（2）编制按从大到小的顺序1211)1(111222-+⋅⋅⋅+--+-=N N S N ，计算N S 的通用程序。

（3）按两种顺序分别计算642101010,,S S S ，并指出有效位数。

（编制程序时用单精度）（4）通过本上机题，你明白了什么？c++代码：#include <iostream.h>#include <fstream.h>void main(){float S1,S2,S=0.0，N,j;fstream myfile; //定义文件myfile.open("chap0.txt",ios::in|ios::out);for (int i=0;i<3;i++){do //读入三组数据{cout<<"请输入N 值:";cin>>N;cout<<endl;myfile<<"当N = "<<N<<" 时"<<endl<<endl;}while (N<2);S=(float) ((1.5-(1.0/N)-(1.0/(N+1)))/2);myfile<<"精确值S=:"<<S<<endl<<endl;S1=0.0;S2=0.0;for (j=2;j<=N;j++) S1=S1+(float) (1.0/(j*j-1.0));for (j=N;j>=2;j--) S2=S2+(float) (1.0/(j*j-1.0));myfile<<"按从大到小的顺序求和结果S1=:"<<S1<<endl<<endl;myfile<<"按从小到大的顺序求和结果S2=:"<<S2<<endl<<endl;}myfile.close();}输出结果：当N = 100 时精确值S=:0.740049按从大到小的顺序求和结果S1=:0.740049 有效位数6按从小到大的顺序求和结果S2=:0.74005 有效位数5当N = 10000 时精确值S=:0.7499按从大到小的顺序求和结果S1=:0.749852 有效位数 4按从小到大的顺序求和结果S2=:0.7499 有效位数 4当N = 1000000 时精确值S=:0.749999按从大到小的顺序求和结果S1=:0.749852 有效位数 4按从小到大的顺序求和结果S2=:0.749999 有效位数 6体会：通过本题，我认识到了舍入误差对计算结果的影响，所以必须选择合适的计算方法，尽量减少这种误差的影响。

数值分析上机题答案

数值分析上机题答案【篇一：数值分析上机试题对应参考答案】么是近似值x* 有效数字？若近似值x*的误差限是某一位的半个单位，该位到x*的第一位非零数字共有n位，就说x*有n位有效数字。

它可表示为2、数值计算应该避免采用不稳定的算法，防止有效数字的损失. 因此，在进行数值运算算法设计过程中主要注意什么？（1）简化计算过程，减少运算次数；（2）避免两个相近的数相减；（3）避免除数的绝对值远小于被除数的绝对值；（4）防止大数“吃掉”小数的现象；（5）使用数值稳定的算法，设法控制误差的传播。

3、写出“n 阶阵a 具有n 个不相等的特征值”的等价条件（至少写3 个）(1)|a|不为零(2)n阶矩阵a的列或行向量组线性无关(3)矩阵a为满秩矩阵(4)n阶矩阵a与n阶可逆矩阵b等价4、迭代法的基本思想是什么？就是用某种极限过程去逐步逼近线性方程组精确解得方法。

其基本思想为：先任取一组近似解初值x0，然后按照某种迭代原则，由x0计算新的近似解x1，以此类推，可计算出x2,x3，…xk，。

,如果｛x｝收敛，则取为原方程组的解。

5、病态线性方程组的主要判断方法有哪些？（1）系数矩阵的某两行（列）几乎近似相关（2）系数矩阵的行列式的值很小（3）用主元消去法解线性方程组时出现小主元（4）近似解x*已使残差向量r=b-ax*的范数很小，但该近似解仍不符合问题要求。

6、lagrange 插值的前提条件是什么？并写出二次lagrange 插值的基函数。

1，j?i?（x)? 前提条件是：l i ,j?0,1,2?，n.?ij0，j?i?二次lagrange 插值的基函数： (x?x)(x?x)12??lx0(xx)(xx) 0?10?2 (x?x)(x?x)02?? lx1(xx)(xx)1?01?2(x?x)(x?x)01?? lx2(x?x)(x?x)20217、什么是数值积分的代数精度？如果某一个求积公式对于次数不超过m的多项式均能准确地成立，但对于m+1次多项式就不准确成立，则称该求积公式具有m次代数精度（或代数精确度）。