数值计算方法大作业
- 格式:doc
- 大小:277.00 KB
- 文档页数:28
数值分析大作业一一、算法设计方案1、求λ1和λ501的值:思路:采用幂法求出按模最大特征值λmax,该值必为λ1或λ501,若λmax小于0,则λmax=λ1;否则λmax=λ501。
再经过原点平移,使用幂法迭代出矩阵A-λmax I的特征值,此时求出的按模最大特征值即为λ1和λ501的另一个值。
2、求λs的值:采用反幂法求出按模最小的特征值λmin即为λs,其中的方程组采用LU分解法进行求解。
3、求与μk最接近的特征值:对矩阵A采用带原点平移的反幂法求解最小特征值,其中平移量为:μk。
4、A的条件数cond(A)=| λmax/λmin|;5、A的行列式的值:先将A进行LU分解,再求U矩阵对角元素的乘积即为A 行列式的值。
二、源程序#include<iostream>#include<iomanip>#include<math.h>#define N 501#define E 1.0e-12 //定义精度常量#define r 2#define s 2using namespace std;double a[N];double cc[5][N];void init();double mifa();double fmifa();int max(int aa,int bb);int min(int aa,int bb);int max_3(int aa,int bb,int cc);void LU();void main(){double a1,a2,d1,d501=0,ds,det=1,miu[39],lamta,cond;int i,k;init();/*************求λ1和λ501********************/a1=mifa();if(a1<0)d1=a1; //若小于0则表示λ1的值elsed501=a1; //若大于0则表示λ501的值for(i=0;i<N;i++)a[i]=a[i]-a1;a2=mifa()+a1;if(a2<0)d1=a2; //若小于0则表示λ1的值elsed501=a2; //若大于0则表示λ501的值cout<<"λ1="<<setiosflags(ios::scientific)<<setprecision(12)<<d1<<"\t";cout<<"λ501="<<setiosflags(ios::scientific)<<setprecision(12)<<d501<<endl;/**************求λs*****************/init();ds=fmifa();cout<<"λs="<<setiosflags(ios::scientific)<<setprecision(12)<<ds<<endl;/**************求与μk最接近的特征值λik**************/cout<<"与μk最接近的特征值λik:"<<endl;for(k=0;k<39;k++){miu[k]=d1+(k+1)*(d501-d1)/40;init();for(i=0;i<N;i++)a[i]=a[i]-miu[k];lamta=fmifa()+miu[k];cout<<"λi"<<k+1<<"\t\t"<<setiosflags(ios::scientific)<<setprecision(12)<<lamta<<en dl;}/**************求A的条件数**************/cout<<"矩阵A的条件式";cond=abs(max(abs(d1),abs(d501))/ds);cout<<"cond="<<setiosflags(ios::scientific)<<setprecision(12)<<cond<<endl;/**************求A的行列式**************/cout<<"矩阵A的行列式";init();LU();for(i=0;i<N;i++){det*=cc[2][i];}cout<<"det="<<setiosflags(ios::scientific)<<setprecision(12)<<det<<endl;system("pause");}/**************初始化函数,给a[N]赋值*************/void init(){int i;for(i=1;i<=501;i++)a[i-1]=(1.64-0.024*i)*sin((double)(0.2*i))-0.64*exp((double)(0.1/i)); }/**************幂法求最大绝对特征值**************/double mifa(){int i,k=0;double u[N],y[N]={0},b=0.16,c=-0.064,Beta_=0,error;for(i=0;i<501;i++)u[i]=1; //令u[N]=1for(k=1;k<2000;k++) //控制最大迭代次数为2000{/***求y(k-1)***/double sum_u=0,gh_sum_u;for(i=0;i<N;i++){sum_u+=u[i]*u[i]; }gh_sum_u=sqrt(sum_u);for(i=0;i<N;i++){y[i]=u[i]/gh_sum_u;}/****求新的uk****/u[0]=a[0]*y[0]+b*y[1]+c*y[2];u[1]=b*y[0]+a[1]*y[1]+b*y[2]+c*y[3]; //前两列和最后两列单独拿出来求中D间的循环求for(i=2;i<N-2;i++){u[i]=c*y[i-2]+b*y[i-1]+a[i]*y[i]+b*y[i+1]+c*y[i+2];}u[N-2]=c*y[N-4]+b*y[N-3]+a[N-2]*y[N-2]+b*y[N-1];u[N-1]=c*y[N-3]+b*y[N-2]+a[N-1]*y[N-1];/***求beta***/double Beta=0;for(i=0;i<N;i++){Beta+=y[i]*u[i];}//cout<<"Beta"<<k<<"="<<Beta<<"\t"; 输出每次迭代的beta /***求误差***/error=abs(Beta-Beta_)/abs(Beta);if(error<=E) //若迭代误差在精度水平内则可以停止迭代{return Beta;} //控制显示位数Beta_=Beta; //第个eta的值都要保存下来,为了与后个值进行误差计算 }if(k==2000){cout<<"error"<<endl;return 0;} //若在最大迭代次数范围内都不能满足精度要求说明不收敛}/**************反幂法求最小绝对特¬征值**************/double fmifa(){int i,k,t;double u[N],y[N]={0},yy[N]={0},b=0.16,c=-0.064,Beta_=0,error;for(i=0;i<501;i++)u[i]=1; //令u[N]=1for(k=1;k<2000;k++){double sum_u=0,gh_sum_u;for(i=0;i<N;i++){sum_u+=u[i]*u[i]; }gh_sum_u=sqrt(sum_u);for(i=0;i<N;i++){y[i]=u[i]/gh_sum_u;yy[i]=y[i]; //用重新赋值,避免求解方程组的时候改变y的值}/****LU分解法解方程组Au=y,求新的***/LU();for(i=2;i<=N;i++){double temp_b=0;for(t=max(1,i-r);t<=i-1;t++)temp_b+=cc[i-t+s][t-1]*yy[t-1];yy[i-1]=yy[i-1]-temp_b;}u[N-1]=yy[N-1]/cc[s][N-1];for(i=N-1;i>=1;i--){double temp_u=0;for(t=i+1;t<=min(i+s,N);t++)temp_u+=cc[i-t+s][t-1]*u[t-1];u[i-1]=(yy[i-1]-temp_u)/cc[s][i-1];}double Beta=0;for(i=0;i<N;i++){Beta+=y[i]*u[i];}error=abs(Beta-Beta_)/abs(Beta);if(error<=E){return (1/Beta);}Beta_=Beta;}if(k==2000){cout<<"error"<<endl;return 0;} }/**************求两数最大值的子程序**************/int max(int aa,int bb){return(aa>bb?aa:bb);}/**************求两数最小值的子程序**************/int min(int aa,int bb){return(aa<bb?aa:bb);}/**************求三数最大值的子程序**************/int max_3(int aa,int bb,int cc){ int tt;if(aa>bb)tt=aa;else tt=bb;if(tt<cc) tt=cc;return(tt);}/**************LU分解**************/void LU(){int i,j,k,t;double b=0.16,c=-0.064;/**赋值压缩后矩阵cc[5][501]**/for(i=2;i<N;i++)cc[0][i]=c;for(i=1;i<N;i++)cc[1][i]=b;for(i=0;i<N;i++)cc[2][i]=a[i];for(i=0;i<N-1;i++)cc[3][i]=b;for(i=0;i<N-2;i++)cc[4][i]=c;for(k=1;k<=N;k++){for(j=k;j<=min(k+s,N);j++){double temp=0;for(t=max_3(1,k-r,j-s);t<=k-1;t++)temp+=cc[k-t+s][t-1]*cc[t-j+s][j-1];cc[k-j+s][j-1]=cc[k-j+s][j-1]-temp;}//if(k<500){for(i=k+1;i<=min(k+r,N);i++){double temp2=0;for(t=max_3(1,i-r,k-s);t<=k-1;t++)temp2+=cc[i-t+s][t-1]*cc[t-k+s][k-1];cc[i-k+s][k-1]=(cc[i-k+s][k-1]-temp2)/cc[s][k-1];}}}}三、程序结果。
数值计算大作业题目一、非线性方程求根1.题目假设人口随时间和当时人口数目成比例连续增长,在此假设下人口在短期内的增长建立数学模型。
(1)如果令()N t 表示在t 时刻的人口数目,β表示固定的人口出生率,则人口数目满足微分方程()()dN t N t dt β=,此方程的解为0()=tN t N e β; (2)如果允许移民移入且速率为恒定的v ,则微分方程变成()()dN t N t vdt β=+, 此方程的解为0()=+(1)t t vN t N e e βββ-;假设某地区初始有1000000人,在第一年有435000人移入,又假设在第一年年底该地区人口数量1564000人,试通过下面的方程确定人口出生率β,精确到410-;且通过这个数值来预测第二年年末的人口数,假设移民速度v 保持不变。
4350001564000=1000000(1)e e βββ+-2.数学原理采用牛顿迭代法,牛顿迭代法的数学原理是,对于方程0)(=x f ,如果)(x f 是线性函数,则它的求根是很容易的,牛顿迭代法实质上是一种线性化方法,其基本思想是将非线性方程0)(=x f 逐步归结为某种线性方程来求解。
设已知方程0)(=x f 有近似根k x (假定0)(≠'x f ),将函数)(x f 在点k x进行泰勒展开,有.))(()()(⋅⋅⋅+-'+≈k k k x x x f x f x f于是方程0)(=x f 可近似地表示为))(()(=-'+k k x x x f x f这是个线性方程,记其根为1k x +,则1k x +的计算公式为)()(1k k k k x f x f x x '-==+,,,2,1,0⋅⋅⋅=k这就是牛顿迭代法,简称牛顿法。
3.程序设计作出函数的图像,大概估计出根的位置fplot('1000*exp(x)+(435*x)*(exp(x)-1)-1564',[0 3]);grid大概估计出初始值x=0.5function [p1,err,k,y]=newton(f,df,p0,delta,max1) % f 是非线性系数 % df 是f 的微商 % p0是初始值% dalta 是给定允许误差 % max1是迭代的最大次数 % p1是牛顿法求得的方程近似解 % err 是p0误差估计 % k 是迭代次数 p0,feval('f',p0) for k=1:max1p1=p0-feval('f',p0)/feval('df',p0); err=abs(p1-p0); p0=p1;p1,err,k,y=feval('f',p1) if(err<delta)|(y==0), break,endp1,err,k,y=feval('f',p1) endfunction y=f(x)y=1000000*exp(x)+435000*(exp(x)-1)/x-1564000; function y=df(x)y=1000000*exp(x)+435000*(exp(x)/x-(exp(x)-1)/x^2);4.结果分析与讨论newton('f','df',1.2,10^(-4),10) 运行后得出结果 p0 =0.5000p1 =0.1679 err =0.3321 k =1 y =9.2415e+004 p1 =0.1031 err =0.0648 k =2 y =2.7701e+003 p1 =0.1010 err =0.0021 k =3 y =2.6953p1 =0.1010 err =2.0129e-006 k =4 y = 2.5576e-006 ans =0.1010运算后的结果为1010.0=β,通过这个数值来预测第二年年末的人口数,0.10100.1010435000f(t)=1000000(1)0.1010t te e +-t=2时候对于f ()2187945.865x =实践表明,当初始值难以确定时,迭代法就不一定收敛了,因此要根据问题实际背景或者二分法先得一个较好的初始值,然后再进行迭代;再者迭代函数选择不合适的话,采用不动点迭代法也有可能出现不收敛的情况;因此我采用的是牛顿法。
题目利用数值计算方法求取基尼系数姓名与学号指导教师年级与专业所在学院一、问题综述:基尼系数(Gini coefficient),是20世纪初意大利学者科拉多·吉尼根据劳伦茨曲线所定义的判断收入分配公平程度的指标。
是比例数值,在0和1之间。
基尼指数(Gini index)是指基尼系数乘100倍作百分比表示。
在民众收入中,如基尼系数最大为“1”,最小等于“0”。
前者表示居民之间的收入分配绝对不平均(即所有收入都集中在一个人手里,其余的国民没有收入),而后者则表示居民之间的收入分配绝对平均,即人与人之间收入绝对平等,但这两种情况只出现在理论上;因此,基尼系数的实际数值只能介于0~1之间,基尼系数越小收入分配越平均,基尼系数越大收入分配越不平均。
设右图中的实际收入分配曲线(红线)和收入分配绝对平等线(绿线)之间的面积为A,和收入分配绝对不平等线(蓝线)之间的面积为B,则表示收入与人口之间的比例的基尼系数为AA+B。
如果A为零,即基尼系数为0,表示收入分配完全平等(红线和绿线重叠);如果B为零,则系数为1,收入分配绝对不平等(红线和蓝线重叠)。
该系数可在0和1之间取任何值。
实际上,一般国家的收入分配,既不是完全平等,也不是完全不平等,而是在两者之间,劳伦茨曲线为一条凸向横轴的曲线。
收入分配越趋向平等,劳伦茨曲线的弧度越小(斜度越倾向45度),基尼系数也越小;反之,收入分配越趋向不平等,劳伦茨曲线的弧度越大,那么基尼系数也越大。
基尼系数的调节需要国家通过财政政策进行国民收入的二次分配,例如对民众的财政公共服务支出和税收等,从而让收入均等化,令基尼系数缩小。
基尼系数由于给出了反映居民之间贫富差异程度的数量界线,可以较客观、直观地反映和监测居民之间的贫富差距,预报、预警和防止居民之间出现贫富两极分化。
因此得到世界各国的广泛认同和普遍采用。
联合国有关组织规定:●若低于0.2表示收入平均;●0.2-0.3表示相对平均;●0.3-0.4表示相对合理;●0.4-0.5表示收入差距大;●0.6以上表示收入差距悬殊。
数值计算方法习题一(2)习题二(6)习题三(15)习题四(29)习题五(37)习题六(62)习题七(70)2009.9,9习题一1.设x >0相对误差为2%,4x 的相对误差。
解:由自变量的误差对函数值引起误差的公式:(())(())'()()()()f x xf x f x x f x f x δδ∆=≈得(1)()f x =11()()*2%1%22x x δδδ≈===;(2)4()f x x =时444()()'()4()4*2%8%x x x x x x δδδ≈===2.设下面各数都是经过四舍五入得到的近似数,即误差不超过最后一位的半个单位,试指出他们各有几位有效数字。
(1)12.1x =;(2)12.10x =;(3)12.100x =。
解:由教材9P 关于1212.m nx a a a bb b =±型数的有效数字的结论,易得上面三个数的有效数字位数分别为:3,4,53.用十进制四位浮点数计算 (1)31.97+2.456+0.1352; (2)31.97+(2.456+0.1352)哪个较精确?解:(1)31.97+2.456+0.1352 ≈21((0.3197100.245610)0.1352)fl fl ⨯+⨯+ =2(0.3443100.1352)fl ⨯+=0.3457210⨯(2)31.97+(2.456+0.1352)21(0.319710(0.245610))fl fl ≈⨯+⨯ = 21(0.3197100.259110)fl ⨯+⨯ =0.3456210⨯易见31.97+2.456+0.1352=0.345612210⨯,故(2)的计算结果较精确。
4.计算正方形面积时,若要求面积的允许相对误差为1%,测量边长所允许的相对误差限为多少?解:设该正方形的边长为x ,面积为2()f x x =,由(())(())'()()()()f x xf x f x x f x f x δδ∆=≈解得(())()()'()f x f x x xf x δδ≈=2(())(())22f x x f x x xδδ==0.5%5.下面计算y 的公式哪个算得准确些?为什么?(1)已知1x <<,(A )11121xy x x-=-++,(B )22(12)(1)x y x x =++; (2)已知1x >>,(A)y =,(B)y = (3)已知1x <<,(A )22sin x y x =,(B )1cos 2xy x-=;(4)(A)9y =(B)y =解:当两个同(异)号相近数相减(加)时,相对误差可能很大,会严重丧失有效数字;当两个数相乘(除)时,大因子(小除数)可能使积(商)的绝对值误差增大许多。
数值计算方法大作业--资料-CAL-FENGHAI-(2020YEAR-YICAI)_JINGBIAN计算方法大作业学生学号: ********学生姓名: ****专业班级: ***********摘要:大作业通过MATLAB在计算方法中的应用实例,探讨了MATLAB在计算方法中的应用方法和技巧,对运用计算机软件完成“计算方法”课程的图形绘制,多项式方程的求解,计算方法分析具有较好的参考价值。
关键字:MATLAB应用迭代法多项式引言在科学研究与工程设计中,经常会遇到数学模型的求解问题,然而在许多情况下,要获得模型问题的准确解是十分困难的,甚至是不可能的。
因此,研究各种数学问题的近似解法非常重要。
数值计算方法又称计算方法或数值计算分析,是一门与计算机应用密切结合的实用性很强的数学课程。
数值计算方法提供的算法具有以下特点:1.面向计算机,根据计算机的特点设计可行的算法。
2.有可靠的理论依据。
3.高效率。
数值计算方法既重视与方法有关的理论,又重视方法的实际运用,而且数值计算方法课程涉及的面较广泛,包括了微积分、线性代数、常微积分方程等数学问题的数值方法。
所以我们只有努力的掌握这几门课程的基本内容,才能学好这门课程。
掌握数值计算方法,包括数组和数组函数,矩阵和矩阵函数的创建与操作,关系与逻辑操作符的运算,多项式计算,数据分析,以及方程与方程组的解法。
掌握Matla图形和3D可视化的技术,围绕数据成图机理,绘图要旨和修饰技法熟悉各种绘图指令和交互操作工具。
包括二维,三维和高维图形绘制,图形的色彩,光源和材质等效果的处理,以及图形句柄操作和动画制作技术。
Matlab数值计算,数值计算功能是Matlab最具代表性的特点,也是最基本、最重要的功能,它是备受欢迎的基石。
Matlab能够成为世界上最优秀的数学软件之一和它出色的数值运算能力是分不开的。
Matlab在数值运算中以数组和矩阵为基础。
数组是Matlab运算中一个重要的数据组织形式。
习 题 一1、下列各数都是经过四舍五入得到的近似值。
试分别指出它们的绝对误差、相对误差和有效数字的位数。
35801=x ,00476.02=x ,33101430.0⨯=x ,24102958-⨯=x ,85000.55=x 。
解:2、已知2031.1=a ,978.0=b 是经过四舍五入后得到的近似值,问b a b a ⨯+,有几位有效数字?解:321110,1022a ab b *-*--≤⨯-≤⨯,而 2.1811, 1.1766a b a b +=⨯=()()3212111101010222a b a b a a b b ****---+-+≤-+-≤⨯+⨯≤⨯故a b +至少具有2位有效数字。
()()32120.978 1.2031110100.006510222ab a b b a a a b b *****----≤-+-≤⨯+⨯=≤⨯ 故a b +至少具有2位有效数字。
3、求二次方程01162=+-x x 的较小正根,取94.763≈,要求有3位有效数字。
解:*1228887.940.06x x x ==-==,*2x 只有一位有效数字。
若改用2180.062715.94x =≈,具有三位有效数字。
4、正方形的边长约cm 100,问测量边长时误差应多大,才能保证面积的误差不超过21cm ?解:正方形的面积函数为2()A x x =(*)2*(*)A A x εε∴=.当*100x =时,若(*)1A ε≤, 则21(*)102x ε-≤⨯故测量中边长误差限不超过0.005cm 时,才能使其面积误差不超过21cm 5、改变下列表达式,使其计算结果比较精确。
(1)1,11>>--+x xx x x ; (2)1),1ln(2>>--x x x ; (3)1,1<<-x e x ; (4)xxsin cos 1-。
解:(1(2) (3)(4)22sin 1cos 2tan sin 22sin cos 22x x xx x x -==习 题 二1、已知314567.032.0sin =,333487.034.0sin =,352274.036.0sin =,用线性插值及抛物插值计算3367.0sin 的值并估计截断误差。
数值计算方法丁丽娟课后习题答案【篇一:北京理工大学数值计算方法大作业数值实验1】)书p14/4分别将区间[?10,10]分为100,200,400等份,利用mesh或surf命令画出二元函数的三维图形。
z=|??|+ ??+?? +??++??【matlab求解】[x,y]=meshgrid(-10:0.1:10);a=exp(-abs(x));b=cos(x+y);c=1./(x.^2+y.^2+1);z=a+b+c;mesh(x,y,z);[x,y]=meshgrid(-10:0.05:10);a=exp(-abs(x));b=cos(x+y);c=1./(x.^2+y.^2+1);z=a+b+c;mesh(x,y,z);[x,y]=meshgrid(-10:0.025:10); a=exp(-abs(x));b=cos(x+y);c=1./(x.^2+y.^2+1);z=a+b+c;mesh(x,y,z);(二)书p7/1.3.2数值计算的稳定性(i)取= ??c语言程序—不稳定解 +=ln1.2,按公式=?? (n=1,2,…) #includestdio.h#includeconio.h#includemath.hvoid main(){float m=log(6.0)-log(5.0),n;int i;i=1;printf(y[0]=%-20f,m); while(i20){n=1/i-5*m;printf(y[%d]=%-20f,i,n);m=n;i++;if (i%3==0) printf(\n); }getch();}(ii) c语言程序—稳定解≈??[ ??+?? +?? ??+??按公式 =??(??)#includestdio.h#includeconio.h#includemath.hvoid main(){float m=(1/105.0+1/126.0)/2,n; k=n,n-1,n-2,…)(【篇二:北京理工大学数值计算方法大作业数值实验4】 p260/1考纽螺线的形状像钟表的发条,也称回旋曲线,它在直角坐标系中的参数方程为= ?????????????????? ?? ??????????= ?????????????? ??曲线关于原点对称,取a=1,参数s的变化范围[-5,5],容许误差限分别是,,和。
第一章作业第一题问题叙述:构造算法并编程序精确计算二次方程的根。
●设a≠0,b2-4ac>0,且有方程ax2+bx+c=0●包括b2≈b2-4ac的情况(a=c=1,b=±1000000.000001)问题分析:对于普通的二次求根公式:x1,2=−b±√b2−4ac2a当b2>>4ac时,分子可能非常小。
由于计算机中的算术运算存在减性抵消的现象,即两个几乎相等的浮点数相减时会引起舍入误差,所以在这种极端条件下用这个公式就会带来很大的误差。
解决方法:1.使用双精度2.使用变换公式x1,2=−2cb±√b2−4ac3.先利用原公式计算较大的根(即分子不会引起减性抵消),再利用公式:x1x2=c a计算较小的根。
问题解决:1.使用双精度:%This program uses double precision to solve the equation of two degree%And make a comparision to the single precisionclear;clc;[a,b,c]=textread('data.txt','%n%n%n'); %read the numbers from data.txtdelta=b*b-4*a*c;x1=(-b+sqrt(delta))/(2*a);x2=(-b-sqrt(delta))/(2*a); %double precisiona=single(a);b=single(b);c=single(c);%use single precisiondelta=single(b*b-4*a*c);x11=single((-b+sqrt(delta))/(2*a));x12=single((-b-sqrt(delta))/(2*a));fid=fopen('out.txt','w');fprintf(fid,'%g %g %g %g',[x1 x2 x11 x12]);fclose(fid);%write the result into the out.txt下面是计算结果:结论:1.在一般情况下,即没有出现b2≈b2−4ac时,无论是单精度还是双精度下均可以得出正确答案。
数值计算第一次大作业实验目的 以Hilbert 矩阵为例,研究处理病态问题可能遇到的困难。
内容 Hilbert 矩阵的定义是,()11/21/31/1/21/31/41/(1)1/31/41/51/(2)1/1/(1)1/(2)1/(21)n i j H h nn n n n n n =⎡⎤⎢⎥+⎢⎥⎢⎥=+⎢⎥⎢⎥⎢⎥++-⎣⎦它是一个对称正定矩阵,而且()n cond H 随着n 的增加迅速增加,其逆矩阵1,()n i j H α-=,这里,2(1)(1)!(1)!(1)[(1)!(1)!]()!()!i j i jn i n j i j i j n i n j α+-+-+-=+----- 1) 画出ln(())~n cond H n 之间的曲线(可以用任何的一种范数)。
你能猜出ln(())~n cond H n 之间有何种关系吗?提出你的猜想并想法验证。
用行范数for n=1:50 for i=1:n for j=1:nA(i,j)=1/(i+j-1);B(i,j)=factorial(n+i-1)*factorial(n+j-1)/((i+j-1)*(factorial(i-1)*factorial(j-1))^2*factorial(n-i )*factorial(n-j));end endresult1=0; for j=1:nresult1=result1+A(1,j); endresult1=log(result1); result2=0; for i=1:n for j=1:nresult2=B(i,j)+result2; endresult(i)=log(result2); endm=max(result);x(n)=result1+m; end plot([1:50],x)对于更大的n 值,由于Hilbert 逆矩阵中的元素过大,溢出,故在此取50以内的n 。
图1 ln(())~n cond H n 关系曲线图猜想ln(())~n cond H n 之间存在线性关系 验证:设ln(()n cond H an b ∞=+ 在以上程序基础上,再添加>>;>> y=x'; >> l=1:40; >> k=l';>> p=polyfit(k,y,1) %一次多项式拟合 p =3.5446 -3.0931% P=polyfit(k,y,2) %二次多项式拟合 p =-0.0008 3.5778 -3.3253 % P=polyfit(k,y,3) %三次多项式拟合0.0000 -0.0033 3.6198 -3.4777% P=polyfit(k,y,4) %四次多项式拟合-0.0000 0.0002 -0.0082 3.6654 -3.5815 % P=polyfit(k,y,5) %五次多项式拟合 p =0.0000 -0.0000 0.0007 -0.0156 3.7107 -3.6542 从上式可以看出,高次项系数相对于一次项和常数项系数要小很多, 所以取ln(() 3.5446 3.0931n cond H n ∞=-2)设D 是n H 的对角线元素开方构成的矩阵。
数值计算方法大作业
嘿,咱今儿来聊聊数值计算方法大作业呀!这可真是个有趣又有点
头疼的事儿呢!
你想想看,数值计算方法就像是一把神奇的钥匙,能打开好多好多
知识的大门。
做数值计算方法大作业的时候,那感觉就好像在探索一
个神秘的宝藏岛,每一步都充满了未知和挑战。
比如说吧,遇到一个复杂的公式,就像是在森林里碰到了一团乱麻,得耐心地一点点解开。
有时候可能会觉得,哎呀,这可咋整呀,咋这
么难呢!但别急呀,咱得静下心来,仔细琢磨。
这不就跟咱平时解一
道特别难的谜题一样嘛,刚开始觉得毫无头绪,可一旦找到那个关键点,嘿,豁然开朗啦!
在做这个大作业的过程中,可千万不能马虎哟!每一个数据都得像
宝贝一样对待,要是不小心弄错了一个,那可能整个结果都跑偏啦!
这就好比盖房子,一块砖没放好,那房子说不定就歪了呀。
而且呀,团队合作也很重要呢!大家一起讨论,一起想办法,那可
比一个人闷头苦干强多啦。
就好像一群小伙伴一起去冒险,每个人都
能发挥自己的长处,互相帮助,多有意思呀!
还有啊,别忘了多检查几遍自己的成果。
这就跟出门前照镜子一样,得看看自己有没有哪里不妥当。
可别嫌麻烦,这可是关乎最后成果好
不好的关键一步呢!
数值计算方法大作业,它既是挑战,也是机会呀!通过完成它,我们能学到好多好多实用的知识和技能,以后遇到类似的问题,咱就可以轻松应对啦,这多棒呀!所以呀,别害怕它,勇敢地去面对,去探索,去享受这个过程吧!咱肯定能把它完成得漂漂亮亮的,让别人都竖起大拇指,你说是不是呢?。
《数值计算方法》在线作业二一、单选题(共 40 道试题,共 100 分。
)1. 欧拉法形式简单,计算方便,但是精度比较低. 正确. 错误正确答案:2. 下列说法中不属于数值方法设计中的可靠性分析的是(). 方法收敛性. 方法稳定性. 方法的计算量. 方法的误差估计正确答案:3. 在插值节点、插值条件相同的情况下,牛顿插值多项式和拉格朗日插值多项式的本质是一样的,只是计算过程不一样. 正确. 错误正确答案:4. 设f(x)=x^4,以-1,0,2,4为节点的三次插值多项式为5x^3-2x^2-8x. 正确. 错误正确答案:5. 迭代法收敛中,若正定,2-也正定,则G-S收敛. 正确. 错误正确答案:6. 设f(0)=0,f(1)=16,f(2)=46,则l1(x)=-x(x-2). 正确. 错误正确答案:7. 用抛物线公式二分前后的两个积分值做线性组合,其结果正好是用科茨公式得到的积分值. 正确. 错误正确答案:8. 区间[,]上的三次样条插值函数是(). 在[,]上2阶可导,节点的函数值已知,子区间上为3次多项式. 在区间[,]上连续的函数. 在区间[,]上每点可微的函数. 在每个子区间上可微的多项式正确答案:9. 解非线性方程f(x)=0的牛顿迭代法在重根附近(). 线性收敛. 三次收敛. 平方收敛. 不收敛正确答案:10. 数值稳定的算法是指舍入误差对计算结果影响不大的算法. 正确. 错误正确答案:11. 如果是正交矩阵,则on2()=(). 0. 1. 2. -1正确答案:12. f(x)=^x在区间[-1,1]上的三次最佳平方逼近多项式为0.9963+0.9979x+0.5367x^2+0.1761x^3. 正确. 错误正确答案:13. 龙贝格积分法是将区间[,]()并进行适当组合而得出的积分近似值的求法. 逐次分半. 回代. 收缩. 拟合正确答案:14. 采用龙格-库塔法求解常微分方程的初值问题时,公式阶数越高,数值解越精确. 正确. 错误正确答案:15.....正确答案:16.....正确答案:17. 预估-校正公式的截断误差为O(h^3). 正确. 错误正确答案:18. 用数值微分公式中求导数值时,步长越小计算就越精确. 正确. 错误正确答案:19. 反幂法用于求矩阵的按模最小的特征值和对应的特征向量,及其求对应于一个给定的近似特征值的特征向量. 正确. 错误正确答案:20. 已知观察值(xi,yi),用最小二乘法求n次拟合多项式Pn(x)时,Pn(x)的次数n可以任意取. 正确. 错误正确答案:21.....正确答案:22.....正确答案:23. 为了使计算y=10+1/(x-1)+2/(x-1)^2-3/(x-1)^3的乘除法运算次数尽量的少,应将表达式改写成y=10+(1+(2-3/(x-1))/(x-1))/(x-1). 正确. 错误正确答案:24. 规格化浮点数系F=(2,4,-1,2)中一共有()个数. 31. 32. 33. 16正确答案:25. 哪种线性方程组可用平方根法求解. 系数矩阵为对称正定. 系数矩阵为正交矩阵. 线性方程组有零解. 线性方程组系数矩阵为单位阵正确答案:26. 迭代法收敛中,若正定,2-非正定,则Joi发散. 正确. 错误正确答案:27.....正确答案:28.....正确答案:29. 下列条件中,不是分段线性插值函数P(x)必须满足的条件为(). P(x)在各节点处可导. P(x)在[,]上连续. P(x)在各子区间上是线性函数. P(xk)=yk正确答案:30.....正确答案:31. 分别改写方程2^x+x-4=0为x=-2^x+4和x=ln(4-x)/ln2的形式,对两者相应迭代公式求所给方程在[1,2]内的实根,下列描述正确的是(). 前者收敛,后者发散. 前者发散,后者收敛. 两者均收敛. 两者均发散正确答案:32. n阶正交矩阵的乘积是()矩阵. 单位. 对称. 实. 正交正确答案:33. 行范数即为无穷范数. 正确. 错误正确答案:34. 设、Q为实空间中矩阵,且有Q^TQ=I,则有||||2=||Q||2. 正确. 错误正确答案:35. Nwton迭代法对于单根是()阶局部收敛的. 一. 二. 三. 四正确答案:36. 使用高斯消去法解线性代数方程组,一般要用选主元的技术. 正确. 错误正确答案:37. 求解常微分方程初值问题的显式欧拉格式具有1阶方法. 正确. 错误正确答案:38. (x-x0)(x-x2)/((x1-x0)(x1-x2))表示在节点x1的二次拉格朗日插值基函数. 正确. 错误正确答案:39. 设数据x1,x2,x3的绝对误差为0.002,那么x1-x2+x3的绝对误差约为0.006. 正确. 错误正确答案:40. 差分方法就是用差商代替微商,降低微分的阶数,以致将微分方程变成代数方程,再用常规方法求解. 正确. 错误正确答案:。
计算方法结课大作业第一题:求数值积分11sin d 3sin xx x x-+⎰,精确到610-。
解:Romberg 方法基本思想: 根据(0)0T (1)0T (0)1T (2)0T (1)1T (0)2T (03)0T (2)1T (1)2T (0)3T 分别算出T 0(0), T 0(1),T 1(0)…….的值,若对角线上的值的差足够小则停止,否则继续计算,直到足够小为止。
T 0(0)=[()()]2b af a f b -+ =11sin(1)sin1[]23sin(1)3sin1+-+-+-+ = sin1sin(1)[]3sin13sin(1)-++-+-其中sin1= 0.841471 sin(-1)= -0.841471=0.841471-0.8414713+0.841471-3+0.841471+=0.438098T 0(1)= 12T 0(0)+1[()]22b a f a b a -+- = 0.469049T 1(0)=(1)(0)00441T T --= 40.469049-0.4380983⨯=0.479366T 0(2)=12T 0(1)2221(21)[()]22n b a n f a b a =--++-∑ =0.476729T1(1)=(2)(1)00441T T--=40.476729-0.4690493⨯=0.479289T2(0)=2(1)(0)112441T T--=160.479289-0.47936615⨯=0.479284T0(3)=12T0(2)4241(21)[()]22nb a nf a b a=--++-∑=0.478645T1(2)=(3)(2)00441T T--=40.478645-0.4767293⨯=0.478645T2(1)=12T0(1)2221(21)[()]22nb a nf a b a=--++-∑=0.438098T3(0)=3(1)(0)223441T T--=640.438098-0.47928463⨯=0.479284化简成缩略图,根据公式可一目了然:T0(0)=0.011568,T0(1)=0.469049, T1(0)= 0.479366T0(2)=0.476729, T1(1)=0.479289 T2(0)=0.479284T0(3)= 0.478645, T1(2)=0.479284, T2(1)=0.479284,T3(0)=0.479284由于T3(0)就已达预定精度,故取I≈T3(0)=0.479284第二题:设实验数据如下建立形如2y ax bx c=++最小二乘拟合。
数值计算方法丁丽娟课后习题答案【篇一:北京理工大学数值计算方法大作业数值实验1】)书p14/4分别将区间[?10,10]分为100,200,400等份,利用mesh或surf命令画出二元函数的三维图形。
z=|??|+ ??+?? +??++??【matlab求解】[x,y]=meshgrid(-10:0.1:10);a=exp(-abs(x));b=cos(x+y);c=1./(x.^2+y.^2+1);z=a+b+c;mesh(x,y,z);[x,y]=meshgrid(-10:0.05:10);a=exp(-abs(x));b=cos(x+y);c=1./(x.^2+y.^2+1);z=a+b+c;mesh(x,y,z);[x,y]=meshgrid(-10:0.025:10); a=exp(-abs(x));b=cos(x+y);c=1./(x.^2+y.^2+1);z=a+b+c;mesh(x,y,z);(二)书p7/1.3.2数值计算的稳定性(i)取= ??c语言程序—不稳定解 +=ln1.2,按公式=?? (n=1,2,…) #includestdio.h#includeconio.h#includemath.hvoid main(){float m=log(6.0)-log(5.0),n;int i;i=1;printf(y[0]=%-20f,m); while(i20){n=1/i-5*m;printf(y[%d]=%-20f,i,n);m=n;i++;if (i%3==0) printf(\n); }getch();}(ii) c语言程序—稳定解≈??[ ??+?? +?? ??+??按公式 =??(??)#includestdio.h#includeconio.h#includemath.hvoid main(){float m=(1/105.0+1/126.0)/2,n; k=n,n-1,n-2,…)(【篇二:北京理工大学数值计算方法大作业数值实验4】 p260/1考纽螺线的形状像钟表的发条,也称回旋曲线,它在直角坐标系中的参数方程为= ?????????????????? ?? ??????????= ?????????????? ??曲线关于原点对称,取a=1,参数s的变化范围[-5,5],容许误差限分别是,,和。
数值计算方法要求: 一、独立完成, 下面五组题目中, 请任选其中一组题目作答, 满分100分;二、 答题步骤:1. 使用A4纸打印学院指定答题纸(答题纸请详见附件);2. 在答题纸上使用黑色水笔....按题目要求手写..作答; 答题纸上全部信息要求手写, 包含中心、 学号、 姓名、 科目、 答题组数等基础信息和答题内容, 请写明题型、 题号;三、 提交方法: 请将作答完成后整页答题纸以图片形式依次粘贴在一个.......Word .... 文档中...上传(只粘贴部分内容图片不给分), 图片请保持正向、 清楚; 1. 上传文件命名为“中心-学号-姓名-科目.doc ” 2. 文件容量大小: 不得超出20MB 。
提醒: 未按要求作答题目........作业..及雷同作业....., .成绩以...0.分记..!题目以下: 第一组:一、 计算题(共56分) 1、 (28分)设有线性方程组b Ax =, 其中⎪⎪⎪⎭⎫⎝⎛=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=582,3015515103531b A(1)求A LU =分解;(2)求方程组解(3)判定矩阵A 正定性2、 (28分)用列主元素消元法求解方程组1311145431221111x x x --⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥-=-⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦二、 叙述题(共44分)1、 (28分)已知方程组Ax b =, 其中1221111,22213A b -⎡⎤⎡⎤⎢⎥⎢⎥==⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦(1)写出该方程组Jacobi 迭代法和Gauss-Seidel 迭代法分量形式;(2)判定(1)中两种方法收敛性, 假如均收敛, 说明哪一个方法收敛愈加快。
2、 (16分)使用高斯消去法解线性代数方程组, 通常为何要用选主元技术?第二组:一、 综合题(共82分)1、 (28分)(1)写出对应三次Lagrange 插值多项式;(2)作均差表, 写出对应三次Newton 插值多项式, 并计算()1.5f 近似值。
目录第一章非线性方程求根 (3)1.1迭代法 (3)1.2牛顿法 (4)1.3弦截法 (5)1.4二分法 (6)第二章插值 (7)2.1线性插值 (7)2.2二次插值 (8)2.3拉格朗日插值 (9)2.4分段线性插值 (10)2.5分段二次插值 (11)第三章数值积分 (13)3.1复化矩形积分法 (13)3.2复化梯形积分法 (14)3.3辛普森积分法 (15)3.4变步长梯形积分法 (16)第四章线性方程组数值法 (17)4.1约当消去法 (17)4.2高斯消去法 (18)4.3三角分解法 (20)4.4雅可比迭代法 (21)4.5高斯—赛德尔迭代法 (23)第五章常积分方程数值法 (25)5.1显示欧拉公式法 (25)5.2欧拉公式预测校正法 (26)5.3改进欧拉公式法 (27)5.4四阶龙格—库塔法 (28)数值计算方法第一章非线性方程求根1.1迭代法程序代码:Private Sub Command1_Click()x0 = Val(InputBox("请输入初始值x0"))ep = Val(InputBox(请输入误差限ep))f = 0While f = 0X1 = (Exp(2 * x0) - x0) / 5If Abs(X1 - x0) < ep ThenPrint X1f = 1Elsex0 = X1End IfWendEnd Sub例:求f(x)=e2x-6x=0在x=0.5附近的根(ep=10-10)1.2牛顿法程序代码:Private Sub Command1_Click()b = Val(InputBox("请输入被开方数x0"))ep = Val(InputBox(请输入误差限ep))f = 0While f = 0X1 = x0 - (x0 ^ 2 - b) / (2 * b)If Abs(X1 - x0) < ep ThenPrint X1f = 1Elsex0 = X1End IfWendEnd Sub例:求56的值。
(ep=10-10)1.3弦截法程序代码:Private Sub Command1_Click()x0 = Val(InputBox("请输入第一个初始值x0"))X1 = Val(InputBox("请输入第二个初始值x1"))ep = Val(InputBox("请输入误差限ep"))f = 0While f = 0X2 = X1 - (X1 ^ 8 - 13) * (X1 - x0) / ((X1 ^ 8 - 13) - (x0 ^ 8 - 13))If Abs(X2 - X1) < ep ThenPrint X2f = 1Elsex0 = X1X1 = X2End IfWendEnd Sub例:求f(x)=x8-13的正根(初始值x1=1,x2=10,ep=10-10)1.4二分法程序代码:Private Sub Command1_Click()a = Val(InputBox("请输入区间端点a"))b = Val(InputBox("请输入区间端点b"))ep = Val(InputBox("请输入误差限ep"))f = 0While f = 0x = (a + b) / 2fx = Exp(-x / 7) * (9 - 2 * x) - 8fa = Exp(-a / 7) * (9 - 2 * a) - 8If fx = 0 Thenf = 1Print "方程的根是", xElseIf fa * fx > 0 Thena = xElseb = xEnd IfIf Abs(b - a) < ep Thenx = (b + a) / 2f = 1Print "方程的根是", xEnd IfEnd IfWendEnd Sub例:求方程f(x)=e-7/x(9-2x)-8在区间[0,1]的实根。
(ep=10-10)第二章插值2.1线性插值程序代码:Private Sub Command1_Click()X0 = Val(InputBox("请输入第一个结点X:"))Y0 = Val(InputBox("请输入第一个结点Y:"))X1 = Val(InputBox("请输入第二个结点X:"))Y1 = Val(InputBox("请输入第二个结点Y:"))f = 0While f = 0x = Val(InputBox("请输入未知点的自变量值X:"))L0 = (x - X1) / (X0 - X1)L1 = (x - X0) / (X1 - X0)y = L0 * Y0 + L1 * Y1Print "x="; x, "y="; yf = Val(InputBox("是否继续(0/1):"))WendEnd Sub例:已知两点(13 , 1)、(49 , 8),求30处的值。
2.2二次插值程序代码:Private Sub Command1_Click()X0 = Val(InputBox("请输入第一个结点X:"))Y0 = Val(InputBox("请输入第一个结点Y:"))X1 = Val(InputBox("请输入第二个结点X:"))Y1 = Val(InputBox("请输入第二个结点Y:"))X2 = Val(InputBox("请输入第三个结点X:"))Y2 = Val(InputBox("请输入第三个结点Y:"))f = 0While f = 0x = Val(InputBox("请输入未知点的自变量值X:"))L0 = (x - X1) * (x - X2) / (X0 - X1) / (X0 - X2)L1 = (x - X0) * (x - X2) / (X1 - X0) / (X1 - X2)L2 = (x - X0) * (x - X1) / (X2 - X0) / (X2 - X1)y = L0 * Y0 + L1 * Y1 + L2 * Y2Print "x="; x, "y="; yf = Val(InputBox("是否继续(0/1):"))WendEnd Sub例:已知三点(81 ,9)、(100 ,10)、(121 ,10),求98处的值。
2.3拉格朗日插值程序代码:Private Sub Command1_Click()Dim x(), y()n = Val(InputBox("请输入插值节点数N"))ReDim x(n), y(n)For i = 0 To nx(i) = Val(InputBox("请输入插值节点x(" + Str(i) + ")"))y(i) = Val(InputBox("请输入插值节点y(" + Str(i) + ")"))Next if = 0While f = 0xx = Val(InputBox("请输入未知点的自变量x:"))Sum = 0For i = 0 To nt = 1For j = 0 To nIf j <> i Thent = t * (xx - x(j)) / (x(i) - x(j))End IfNext jSum = Sum + t * y(i)Next iPrint "x="; xx, "y="; Sumf = Val(InputBox("是否继续(0/1)"))WendEnd Sub例:已知四点(100 ,10)、(81 ,9)、(64 ,8)、(49 ,7),求87 处的值。
2.4分段线性插值程序代码:Private Sub Command1_Click()Dim x(), y()n = Val(InputBox("请输入插值节点数N"))ReDim x(n), y(n)For i = 0 To nx(i) = Val(InputBox("请输入插值节点x(" + Str(i) + ")")) y(i) = Val(InputBox("请输入插值节点y(" + Str(i) + ")")) Next if = 0While f = 0xx = Val(InputBox("请输入未知点的自变量x:"))L = 0j = 1While L = 0If xx < x(j) Thenk = j + 1L = 1Elsej = j + 1If j > n - 1 Thenk = n - 1L = 1End IfEnd IfWendl0 = (xx - x(k)) / (x(k - 1) - x(k))l1 = (xx - x(k - 1)) / (x(k) - x(k - 1))yy = l0 * y(k - 1) + l1 * y(k)Print "x="; xx, "y="; yyf = Val(InputBox("是否继续(0/1)"))WendEnd Sub例:已知三点(361 , 19)、(324 ,18)、(289 ,17),N=2,求300处的值。
2.5分段二次插值程序代码:Private Sub Command1_Click()Dim x(), y()n = Val(InputBox("请输入插值节点数N"))ReDim x(n), y(n)For i = 0 To nx(i) = Val(InputBox("请输入插值节点x(" + Str(i) + ")"))y(i) = Val(InputBox("请输入插值节点y(" + Str(i) + ")"))Next if = 0While f = 0xx = Val(InputBox("请输入未知点的自变量x:"))If x0 < x(1) Thenk = 1f = 1End Ifi = 2Do While f = 0 And i >= n - 1If x0 < x(i) ThenIf x0 - x(i - 1) < x(i) - x0 Thenk = i - 1f = 1Elsek = if = 1End IfElsei = i + 1End IfLoopIf f = 0 Thenk = n - 1End Ifl1 = (xx - x(k + 1)) * (xx - x(k)) / ((x(k - 1) - x(k + 1)) * (x(k - 1) - x(k)))l2 = (xx - x(k + 1)) * (xx - x(k - 1)) / ((x(k) - x(k + 1)) * (x(k) - x(k - 1)))l3 = (xx - x(k)) * (xx - x(k - 1)) / ((x(k + 1) - x(k)) * (x(k + 1) - x(k - 1)))yy = l1 * y(k - 1) + l2 * y(k) + l3 * y(k + 1)Print "x="; xx, "y="; yyf = Val(InputBox("是否继续(0/1)"))WendEnd Sub例:已知三点(225 , 15)、(196 ,14)、(169 ,13),求180处的值。