数值分析上机题Matlab(东南大学)3
MATLAB 上机 习题及答案

15、今有多项式P1(x)=x4-2x+1,P2(x)=x2+4x-0.5,要求先求得P(x)=P1(x)+P2(x),然后计算xi=0.2*i各点上的P(xi)(i=0,1,2,…,5)值。
p1=[1.0 0.0 0.0 -2.0 1.0];>> p2=[0.0 0.0 1.0 4.0 -0.5];>> p1x=poly2sym(p1);p2x=poly2sym(p2);>> p=p1x+p2xp =x^4+2*x+1/2+x^2>> x=0:5;>> x.^4+2*x+1/2+x.^2ans =0.5000 4.5000 24.5000 96.5000 280.5000 660.50001、试个MATLAB的工作空间中建立以下2个矩阵:A=[1 2]1234B⎡⎤=⎢⎥⎣⎦,求出矩阵A和B的乘积,并将结果赋给变量C。
>> A=[1 2]A =1 2>> B=[1 23 4]B =1 23 4>> C=A*BC =7 102、利用MATLAB提供的帮助信息,了解inv命令的调用格式,并作简要说明。
help invINV Matrix inverse.INV(X) is the inverse of the square matrix X.A warning message is printed if X is badly scaled ornearly singular.See also SLASH, PINV, COND, CONDEST, LSQNONNEG, LSCOV. Overloaded methodshelp gf/inv.mhelp zpk/inv.mhelp tf/inv.mhelp ss/inv.mhelp lti/inv.mhelp frd/inv.mhelp sym/inv.mhelp idmodel/inv.m3、使用help命令查询函数plot的功能以及调用方法,然后利用plot命令绘制函数y=sin(x)的图形,其中0xπ≤≤。
数值分析matlab版第三章

实验报告三题目:曲线拟合的最小二乘法摘要:在实验中往往要从一组实验数据(xi,yi)(i=0,1,2,3,...,m)中寻找自变量x与因变量y之间的函数关系y=F(x)。
由于观测数据往往不准确,因为要求y=F (x)经过所有点(xi,yi),而只要求在给定点xi上误差按某种标准最小。
数学原理:所谓的曲线拟合的最小二乘法就是使拟合多项式在各节点处的偏差P(xi)-yi的平方和2)^)((pn1iyii x-∑=达到最小。
习题22.function c =myployfit(x,y,m) n =length(x);b =zeros(1,m+1);f =zeros(n,m+1);for k=1:m+1f(:,k)=x'.^(k-1);enda=f'*f;b=f'*y';c=a\b;c=flipud(c);>> x=19:6:44;>> y=[19.0,32.3,49.0,73.3,97.8];>>m=1;>> C=mypolyfit(x,y,m)C =3.309999999999997-48.329999999999913>> plot(x,y,'ok','LineWidth',2);>> grid on;>> title('最小二乘一次拟合效果图'); >> x0=0:2:44;>> y0=C(1).*x0+C(2);>> hold on,plot(x0,y0,'-r');051015202530354045-50050100>> x=19:6:44;>> y=[19.0,32.3,49.0,73.3,97.8];>> m=2;>> C=mypolyfit(x,y,m)C =0.059523809523808-0.3804761904760984.586666666665309>> figure;>> plot(x,y,'ok','LineWidth',2);>> grid on;>> title('最小二乘二次拟合效果图');>> x0=0:2:44;>>>> y0=C(1).*x0.^2+C(2).*x0+C(3);>> hold on ,plot(x0,y0,'-r');051015202530354045020406080100120习题23function c =myployfit(t,s,m)n =length(t);b =zeros(1,m+1);f =zeros(n,m+1);for k=1:m+1f(:,k)=t'.^(k-1);enda=f'*f;b=f'*s';c=a\b;c=flipud(c);>> t=[0,0.9,1.9,3.0,3.9,5.0];>> s=[0,10,30,50,80,110];>> m=1;>> C=mypolyfit(t,s,m)C =22.253760999148451-7.855047781247040>> plot(t,s,'ok','LineWidth',2);>> grid on;>> title('最小二乘一次拟合效果图');>> t0=0:5;>> s0=C(1).*t0+C(2);>> hold on ,plot(t0,s0,'-r');00.51 1.52 2.53 3.54 4.55-20020406080100120最小二乘一次拟合效果图>> m=2;>> C=mypolyfit(t,s,m)C =2.24880969100440011.081396145687567-0.583364511695529>> figure;>> plot(t,s,'ok','LineWidth',2);>> grid on;>> title('最小二乘二次拟合效果图');>> t0=0:5;>> y0=C(1).*x0.^2+C(2).*x0+C(3);>> hold on ,plot(x0,y0,'-r');00.51 1.52 2.53 3.54 4.55-20020406080100120最小二乘二次拟合效果图function c =myployfit(t,y,m)n =length(t);b =zeros(1,m+1);f =zeros(n,m+1);for k=1:m+1f(:,k)=t'.^(k-1);enda=f'*f;b=f'*y';c=a\b;c=flipud(c);习题24:>> t=0:5:55;>> y=[0,1.27,2.16,2.86,3.44,3.87,4.15,4.37,4.51,4.58,4.62,4.64];>> m=2;>> C=mypolyfit(t,y,m)C =-0.0022306693306690.1992252747252740.245302197802203>> plot(t,y,'ok','LineWidth',2);>> plot(t,y,'ok','LineWidth',2);>> grid on;>> title('最小二乘一次拟合效果图');>> t0=0:5:55;>> y0=C(1).*t0.^2+C(2).*t0+C(3);>> hold on ,plot(t0,y0,'-r');总结:曲线拟合的最小二乘法是计算机上有效的和节约时间的算法。
matlab上机习题详细讲解_试题答案

P第一次实验答案要求以0.01秒为间隔,求出y的151个点, 并求出其导数的值和曲线。
clcclearx=0:0.01:1.5;y=sqrt(3)/2*exp(-4*x).*si n(4*sqrt(3)*x+pi/3) y1=diff(y)subplot(2,1,1)Plot(x,y)subplot(2,1,2) plot(x(1:150),y1)2绘制极坐标系下曲线(a,b,n自定数据)— a cos b n vclccleara=10;b=pi/2;n=5;theta=0:pi/100:2*pi; rho=a*cos(b+n*theta);polar(theta,rho) z2=X.*2-Y.*3;xlabel( 'x')ylabel( 'y')zlabel( 'z') surf(X,Y,z1)hold on surf(X ,Y, z2)k=fi nd(abs(z1-z2)<0.5);x1=X(k)y1=Y(k) z3=x1.A2-2*y1.A2 hold onplot3(x1,y1,z3, '*')4、设y cos x 0.53sin x(1 x2)把x=0~2 n间分为101点,画出以x为横坐标,y 为纵坐标的曲线,要求有图形标注clcclearx=-2*pi:0.1: 2*pi;y=cos(x).*(0.5+si n(x)*3./(1+x.A2));plot(x,y, 'b*-');title('绘图’);xlabel( 'x 坐标');ylabel( 'y 坐标');legend('原函数')gtext( 'y=cos(x)(0.5+3*sin(x)/(1+xA2))' )3.列出求下列空间曲面交线的程序乙=x2 _2y2z2 = 2x _ 3yclcclearx=[-5:0.5:5];[X,Y]=meshgrid(x);z1=X.A2-2*Y.A2; 5、求下列联立方程的解3x 4y - 7z - 12w 二45x - 7y 4z 2w - -3x 8z - 5w = 9-6x 5y - 2z 10w = -8clccleara=[3,4,-7,-12;5,-7,4,2;1,0,8,-5;-6,5,-2,10]; b=[4,-3,9,-8];第二次试验答案1、编制m文件,等待键盘输入,输入密码20120520 ,密码正确,显示输入密码正确,程序结束;否则提示,重新输入。
东南大学数值分析上机题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]);}。
数值分析上机实习题及答案.docx

方詡文金兴:爭[数值分析]2017-2018第二学期上机实习题1:编程计算亍丄,其中C= 4. 4942x10307,给出并观察计算结心C"果,若有问题,分析之。
解:mat lab 编程如下:E) funct ion diy i ti formatlong g;n 二input C 输入ii 值= c= 4.4942E307; sum 二0; s 二 0;E3 for i = l:n s = l/ (c*i);>> diyiti 输入n 值:10 104.6356e-308 >> diyiti输入ri 值:1001004.6356e-308 >> diyiti 输入n 值:1000 10004.6356e-308 >> diyiti揄入n 值* 1000001000004・ 6356e-308 >> diyiti输入n 值;1000000001000000004.6356e-308图二:运行结果Mat lab 中,forma t long g 对双精度,显示15位定点或浮点格式,由上图 可知,当输入较小的n 值5分别取10, 100, 1000, 100000, 100000000)的时候, 结果后面的指数中总是含有e-308,这和题目中的C 值很相似,我认为是由于分 母中的C 值相对于n 值过大,出现了 “大数吃小数”的现彖,这是不符合算法原 则的。
2:利用牛顿法求方程-1^ = 2于区间241的根,考虑不同初值下牛顿法的收敛情况。
解:牛顿法公式为:利用mat lab 编程function di2ti21 3i=l ;2 2.85208156699784 xO 二input ('输入初值x0:‘ );A 二[i x0];3 2.55030468822809 t=x0+ (x0-log (xO) -2) /(1-1/xO) ; %迭代函数4 1.91547247100476 三 while (abs (t _x0)>0.01)i=i+l; 5 0.37867158538991 xO 二 t; 6 0.774964959780275 A = [A;i xO];t =x0+(x0-log(xO)-2)/(1-1/xO): 7 4.11574081641933 cnd| 8 5.04162436446126 disp (A);96.81782826645596当输入初值二3的时候并不能收敛。
东南大学数值分析上机实验题(下)

数值分析上机报告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)本次上机作业让我知道了在遇到复杂问题中,无法给出解析式的情况下,要学会灵活使用各种微分数值解法,并且能计算出不同方法的精度大小。
matlab上机习题答案
matlab上机习题答案Matlab上机习题答案在现代科学和工程领域中,计算机编程和数值计算已经成为必不可少的技能。
而Matlab作为一种强大的数值计算软件,被广泛应用于各种领域。
为了帮助学习者更好地掌握Matlab的基本操作和数值计算方法,老师们经常会布置一些上机习题,让学生通过实际操作来加深对Matlab的理解。
下面我们来看一些常见的Matlab上机习题答案:1. 编写一个Matlab程序,计算并输出1到100之间所有奇数的和。
答案:```matlabsum = 0;for i = 1:2:100sum = sum + i;enddisp(sum);```2. 编写一个Matlab程序,计算并输出斐波那契数列的前20个数字。
答案:```matlabfib = zeros(1,20);fib(1) = 1;fib(2) = 1;for i = 3:20fib(i) = fib(i-1) + fib(i-2);enddisp(fib);```3. 编写一个Matlab程序,求解一元二次方程ax^2 + bx + c = 0的根。
答案:```matlaba = 1;b = -3;c = 2;delta = b^2 - 4*a*c;if delta < 0disp('无实根');elseif delta == 0x = -b / (2*a);disp(x);elsex1 = (-b + sqrt(delta)) / (2*a);x2 = (-b - sqrt(delta)) / (2*a);disp(x1);disp(x2);end```通过以上几个例子,我们可以看到,Matlab的语法简洁明了,功能强大。
通过编写程序来解决实际问题,不仅加深了对Matlab的理解,也提高了计算机编程和数值计算的能力。
希望大家在学习Matlab的过程中能够多多练习,不断提高自己的编程水平。
Matlab上机题及答案
1 一个三位整数各位数字的立方和等于该数本身则称该数为水仙花数。
输出全部水仙花数。
for m=100:999m1=fix(m/100); %求m的百位数字m2=rem(fix(m/10),10); %求m的十位数字m3=rem(m,10); %求m的个位数字if m==m1*m1*m1+m2*m2*m2+m3*m3*m3disp(m)endend2.从键盘输入若干个数,当输入0时结束输入,求这些数的平均值和它们之和。
sum=0;n=0;val=input('Enter a number (end in 0):');while (val~=0)sum=sum+val;n=n+1;val=input('Enter a number (end in 0):');endif (n > 0)summean=sum/nend3. 若一个数等于它的各个真因子之和,则称该数为完数,如6=1+2+3,所以6是完数。
求[1,500]之间的全部完数。
for m=1:500s=0;for k=1:m/2if rem(m,k)==0s=s+k;endendif m==sdisp(m);endend4. 从键盘上输入数字星期,在屏幕上显示对应英文星期的单词。
function weekn=input('input the number:');if isempty(n)errror('please input !!')endif n>7|n<1error('n between 1 and 7')endswitch ncase 1disp('Monday')case 2disp('Tuesday')case 3disp('Wednesday')case 4disp('Thursday')case 5disp('Friday')case 6disp('Saturday')case 7disp('Sunday')end5. 某公司销售电脑打印机的价格方案如下:()如果顾客只买一台打印机,则一台的基本价格为$150。
数值分析matlab上机作业报告
一、给定向量x ≠0,计算初等反射阵H k 。
1.程序功能:给定向量x ≠0,计算初等反射阵H k 。
2.基本原理: 若()xx x R x ∈=,,, 的分量不全为零,则由12112212122()x (,,,)1()22n T T sign x e x x x x σσσρσσρ-⎧=⎪=+=+⎪⎪⎪==+⎨⎪⎪=-=-⎪⎪⎩u x u uu H I I uuu 确定的镜面反射阵H 使得y e Hx =-=σ;当(1)k n ≤<时,由21/2k ()T 1()()()k 1()()()(())(0,,0,,,,)1()()=()2()nk i i kk nk k k n k T k k Tk k k kk k T k k sign x x x x x x σσρσσσρ=+-⎧=⎪⎪⎪=+∈⎨⎪==+⎪⎪=-⎩∑u R u u u H I u u 有T 121(,,,,,0,,0)n k k k x x x σ-=-∈H x R 算法:(1)输入x ,若x 为零向量,则报错 (2)将x 规范化,{}x x x M ,,,max =如果M =0,则报错同时转出停机 否则n i M x x i i ,,2,1, =←(3)计算2x =σ,如果0<1x ,则σσ-= (4))(1x +=σσρ (5)计算1,(1)x σ==+u x u (6)1Tρ-=-H I uu (7)(M ,0,,0)σ=-y(8)按要求输出,结束3.变量说明:x -输入的n维向量;n -n维向量x的维数;M -M是向量x的无穷范数,即x中绝对值最大的一项的绝对值;p -Householder初等变换阵的系数ρ;u -Householder初等变换阵的向量Us -向量x的二范数;x -输入的n维向量;n -n维向量x的维数;p -Householder初等变换阵的系数ρ;u -Householder初等变换阵的向量Uk -数k,H*x=y,使得y的第k+1项到最后项全为零;4.程序代码:(1)function [p,u]=holder2(x)%HOLDER2 给定向量x≠0,计算Householder初等变换阵的p,u%程序功能:函数holder2给定向量x≠0,计算Householder初等变换阵的p,u;%输入:n维向量x;%输出:[p,u]。
数值分析上机题参考答案.docx
如有帮助欢迎下载支持数值分析上机题姓名:陈作添学号: 040816习题 120.(上机题)舍入误差与有效数N11 3 1 1设S N,其精确值为 。
22 2 NN 1j 2j1(1)编制按从大到小的顺序111 ,计算 S 的通用程序。
S N1 321N 21 N22(2)编制按从小到大的顺序111,计算 S 的通用程序。
S N1(N 1)2 122 1NN 2(3)按两种顺序分别计算S 102 , S 104 , S 106 ,并指出有效位数。
(编制程序时用单精度)(4)通过本上机题,你明白了什么?按从大到小的顺序计算 S N 的通用程序为: 按从小到大的顺序计算 S N 的通用程序为:#include<iostream.h> #include<iostream.h> float sum(float N) float sum(float N) {{float j,s,sum=0; float j,s,sum=0; for(j=2;j<=N;j++) for(j=N;j>=2;j--) {{s=1/(j*j-1); s=1/(j*j-1); sum+=s;sum+=s;}}return sum;return sum;}}从大到小的顺序的值从小到大的顺序的值精确值有效位数从大到小从小到大0.7400490.740050.74004965 S 1020.7498520.74990.74994 4S 1040.7498520.7499990.74999936S 106通过本上机题, 看出按两种不同的顺序计算的结果是不相同的,按从大到小的顺序计算的值与精确值有较大的误差, 而按从小到大的顺序计算的值与精确值吻合。
从大到小的顺序计算得到的结果的有效位数少。
计算机在进行数值计算时会出现“大数吃小数”的现象,导致计算结果的精度有所降低,我们在计算机中进行同号数的加法时,采用绝对值较小者先加的算法,其结果的相对误差较小。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
数值分析上机题2 2010/9/29 一、题目 P.113-114, 35.(上机题)列主元Gauss消去法、逐次超松弛迭代法
二、程序 35 - (1) 解方程根的通用程序如下(Matlab): function x=getGausssResult(A,b) n = size(A,1); exA = [A,b]; % 若方程组无唯一解,返回 0 if rank(A) < n x = 0; return; end % Gauss列主元消去法 for i = 1:n maxx = 0; flag = 0; % 找到第 i 列中值最大的行 for k = i:n if abs(exA(k,i)) > maxx maxx = abs(exA(k,i)); flag = k; end end % 将第 i 列中值最大的行与第 i 行交换 if flag ~= i for k = i:n+1 tmp = exA(i,k); exA(i,k) = exA(flag,k); exA(flag,k) = tmp; end end for j = i+1:n if exA(j,i) ~= 0 l = exA(j,i)/exA(i,i); for k = i:n+1 exA(j,k) = exA(j,k) - l*exA(i,k); end end end end x = zeros(n,1); for i = n:-1:1 tmp = 0; for j = i+1:n tmp = tmp + exA(i,j)*x(j); end x(i) = (exA(i,n+1)-tmp) / exA(i,i); end end
35 - (2) 求解: function exp3_35 global R b R = [ 31 -13 0 0 0 -10 0 0 0;... -13 35 -9 0 -11 0 0 0 0;... 0 -9 31 -10 0 0 0 0 0;... 0 0 -10 79 -30 0 0 0 -9;... 0 0 0 -30 57 -7 0 -5 0;... 0 0 0 0 -7 47 -30 0 0;... 0 0 0 0 0 -30 41 0 0;... 0 0 0 0 -5 0 0 27 -2;... 0 0 0 -9 0 0 0 -2 29]; b = [ -15 27 -23 0 -20 12 -7 7 10]'; x = getGausssResult(R,b); fprintf('========================================'); fprintf('========================================'); fprintf('========================================\n'); fprintf('习题 3_35\n'); fprintf('========================================'); fprintf('========================================'); fprintf('========================================\n'); for i = 1:9 fprintf('\t\tx%d\t',i); end fprintf('\n----------------------------------------'); fprintf('----------------------------------------'); fprintf('----------------------------------------\n'); for i = 1:9 fprintf('\t% 10.5f',x(i)); end fprintf('\n========================================'); fprintf('========================================'); fprintf('========================================\n'); end
36 - (1) 通用程序如下 function [x,n]=getSORResult(x,omega,epsilon) global b D UTilde LTilde SOmega = inv(D+omega*LTilde)*((1-omega)*D-omega*UTilde); fOmega = omega*inv(D+omega*LTilde)*b; xlast = x; x=SOmega*xlast+fOmega; n = 1; while norm(x-xlast)>epsilon xlast = x; x=SOmega*xlast+fOmega; n = n + 1; end end 36 - (2) 程序运行结果见下。三、结果 ======================================================================================================================== 习题 3_35 ======================================================================================================================== x1 x2 x3 x4 x5 x6 x7 x8 x9 ------------------------------------------------------------------------------------------------------------------------ -0.28923 0.34544 -0.71281 -0.22061 -0.43040 0.15431 -0.05782 0.20105 0.29023 ========================================================================================================================
======================================================================================================================================================================== 习题 3_36 ======================================================================================================================================================================== Omega n x1 x2 x3 x4 x5 x6 x7 x8 x9 ------------------------------------------------------------------------------------------------------------------------------------------------------------------------ 0.02 1420 -0.28883 0.34580 -0.71259 -0.22046 -0.43021 0.15472 -0.05740 0.20111 0.29030 0.04 764 -0.28904 0.34562 -0.71271 -0.22054 -0.43031 0.15451 -0.05762 0.20108 0.29026 0.06 528 -0.28910 0.34555 -0.71274 -0.22056 -0.43034 0.15444 -0.05769 0.20107 0.29025 0.08 405 -0.28914 0.34552 -0.71276 -0.22057 -0.43036 0.15441 -0.05772 0.20107 0.29025 0.10 328 -0.28916 0.34551 -0.71277 -0.22058 -0.43037 0.15439 -0.05774 0.20106 0.29024 0.12 276 -0.28917 0.34549 -0.71278 -0.22059 -0.43037 0.15437 -0.05776 0.20106 0.29024 0.14 238 -0.28918 0.34548 -0.71278 -0.22059 -0.43038 0.15436 -0.05777 0.20106 0.29024 0.16 209 -0.28919 0.34548 -0.71279 -0.22059 -0.43038 0.15435 -0.05778 0.20106 0.29024 0.18 186 -0.28919 0.34547 -0.71279 -0.22059 -0.43038 0.15435 -0.05778 0.20106 0.29024 0.20 168 -0.28920 0.34547 -0.71279 -0.22060 -0.43038 0.15434 -0.05779 0.20106 0.29023