数值计算方法大作业
- 格式:doc
- 大小:204.00 KB
- 文档页数:30
数值计算方法练习题习题一1. 下列各数都是经过四舍五入得到的近似数,试指出它们有几位有效数字以及它们的绝对误差限、相对误差限。
(1);(2);(3);(4);(5);(6);(7);2. 为使下列各数的近似值的相对误差限不超过,问各近似值分别应取几位有效数字?3. 设均为第1题所给数据,估计下列各近似数的误差限。
(1);(2);(3)4. 计算,取,利用下列等价表达式计算,哪一个的结果最好?为什么?(1);(2);(3)(4)5. 序列满足递推关系式若(三位有效数字),计算时误差有多大?这个计算过程稳定吗?6. 求方程的两个根,使其至少具有四位有效数字(要求利用。
7. 利用等式变换使下列表达式的计算结果比较精确。
(1);(2)(3);(4)8. 设,求证:(1)(2)利用(1)中的公式正向递推计算时误差增大;反向递推时误差函数减小。
9.设x>0,x*的相对误差为δ,求f(x)=ln x的误差限。
10.下列各数都是经过四舍五入得到的近似值,试指出它们有几位有效数字,并给出其误差限与相对误差限。
11.下列公式如何才比较准确?(1)(2)12.近似数x*=0.0310,是位有数数字。
13.计算取,利用式计算误差最小。
四个选项:习题二1. 已知,求的二次值多项式。
2. 令求的一次插值多项式,并估计插值误差。
3. 给出函数的数表,分别用线性插值与二次插值求的近似值,并估计截断误差。
0.4 0.5 0.6 0.7 0.80.38942 0.47943 0.56464 0.64422 0.717364. 设,试利用拉格朗日余项定理写出以为节点的三次插值多项式。
5. 已知,求及的值。
6. 根据如下函数值表求四次牛顿插值多项式,并用其计算和的近似值。
X 1.615 1.634 1.702 1.828 1.921F (x) 2.41450 2.46459 2.65271 3.03035 3.340667. 已知函数的如下函数值表,解答下列问题(1)试列出相应的差分表;(2)分别写出牛顿向前插值公式和牛顿向后插值公式。
数值分析大作业一一、算法设计方案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分,共17分)1、如果用二分法求方程043=-+x x 在区间]2,1[内的根精确到三位小数,需对分( )次。
2、迭代格式)2(21-+=+k k k x x x α局部收敛的充分条件是α取值在( )。
3、已知⎪⎩⎪⎨⎧≤≤+-+-+-≤≤=31)1()1()1(2110)(233x c x b x a x x x x S 是三次样条函数,则a =( ),b =( ),c =( )。
4、)(,),(),(10x l x l x l n 是以整数点n x x x ,,,10 为节点的Lagrange 插值基函数,则∑==nk kx l0)(( ),∑==nk k jk x lx 0)(( ),当2≥n 时=++∑=)()3(204x l x xk k n k k( )。
5、设1326)(247+++=x x x x f 和节点,,2,1,0,2/ ==k k x k 则=],,,[10n x x x f 和=∆07f。
6、5个节点的牛顿-柯特斯求积公式的代数精度为 ,5个节点的求积公式最高代数精度为 。
7、{}∞=0)(k kx ϕ是区间]1,0[上权函数x x =)(ρ的最高项系数为1的正交多项式族,其中1)(0=x ϕ,则⎰=14)(dx x x ϕ 。
8、给定方程组⎩⎨⎧=+-=-221121b x ax b ax x ,a 为实数,当a 满足 ,且20<<ω时,SOR 迭代法收敛。
9、解初值问题00(,)()y f x y y x y '=⎧⎨=⎩的改进欧拉法⎪⎩⎪⎨⎧++=+=++++)],(),([2),(]0[111]0[1n n n n n n n n n n y x f y x f h y y y x hf y y 是阶方法。
10、设⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=11001a a a a A ,当∈a ( )时,必有分解式T LL A =,其中L 为下三角阵,当其对角线元素)3,2,1(=i l ii 满足( )条件时,这种分解是唯一的。
1.解方程组:x1+x2+x3=6;x1+3x2-2x3=1;2x1-2x2+x3=1.x1,x2,x3分别为()A.1,1,1B.1,1,2C.1,2,3D.3,2,1【参考答案】: C2.若A是n*n阶非奇异阵,则必存在单位下三角阵L和上三角阵U,使A=LU唯一成立A.正确B.错误【参考答案】: B3.设Ax=b,准确解为X*,某一近似解为X,用()来判断误差A.||AX-b||B.||X-X*||C.bD.||b-AX||【参考答案】: B4.已知多项式P(x),过点(0,0)(2,8)(4,64)(11,1331)(15,3375),它的三阶差商为常数1,一阶二阶差商均不是0,那么P(x)是()A.二次多项式B.不超过二次的多项式C.三次多项式D.四次多项式【参考答案】: C5.三次样条函数是在各个子区间上的()次多项式A.2B.3C.4D.5【参考答案】: B6.若误差限为0.5×10^(-5),那么近似数0.003400有5位有效数字A.正确B.错误【参考答案】: B7.A.AB.BC.CD.D【参考答案】: C8.区间[a,b]上的三次样条函数是一个次数不超过三次的多项式A.正确B.错误【参考答案】: B9.数值求积公式中的Simpson公式的代数精度为A.0B.1C.2D.3【参考答案】: D10.如果插值结点相同,在满足相同插值条件下所有的插值多项式是等价的。
A.正确B.错误【参考答案】: A11.A.AB.BC.CD.D【参考答案】: A12.下列说法错误的是()A.如果一个近似数的每一位都是有效数字,则称该近似数为有效数B.凡是经“四舍五入”得到的近似数都是有效数C.数值方法的稳定性是指初始数据的扰动对计算结果的影响D.病态问题是由数学问题本身的性质决定的,与数值方法有关【参考答案】: B13.A.AB.BC.CD.D【参考答案】: D14.双点割线法是用过两点的割线与x轴交点的横坐标作为X0的新近似值A.正确B.错误【参考答案】: A15.设求方程f(x)=0的根的牛顿法收敛,则它具有()收敛A.超线性B.平方C.线性D.三次【参考答案】: C16.经过A(0,1),B(1,2),C(2,3)的插值多项式P(x)=()A.xB.x1C.2x1D.x^21【参考答案】: D17.A.AB.BC.CD.D【参考答案】: A18.当x=1,-1,2时,f(x)=0,-3,4,则f(x)的二次插值多项式为5x^2/6+2x/3-7/3A.正确B.错误【参考答案】: A19.拉格朗日插值基函数在节点上的取值是(0或1?)A.正确B.错误【参考答案】: A20.x1=1,f(x1)=4;x2=2,f(x2)=1;x3=4,f(x3)=6.则二阶差商f(x1,x2,x3)=()A.5/2B.-3C.11/6D.2/5【参考答案】: C21.含有n+1个节点的插值型数值积分公式的代数精度至少为()A.1B.nC.n1D.2【参考答案】: B22.若线性方程组的系数矩阵A的条件数Cond(A)=1,则称此方程为良态方程组A.正确B.错误【参考答案】: B23.数值计算中主要研究的误差有()A.模型误差和观测误差B.截断误差和方法误差C.相对误差和绝对误差 D.模型误差和方法误差【参考答案】: C24.下列求积公式中用到外推技术的是()A.梯形公式B.复合抛物线公式C.龙贝格公式D.高斯型求积公式【参考答案】: B25.下列说法不正确的是()A.二分法不能用于求函数f(x)=0的复根B.方程求根的迭代解法的迭代函数为?f(x),则迭代收敛的充分条件是?f(x)1C.用高斯消元法求解线性方程组AX=B时,在没有舍入误差的情况下得到的都是精确解D.如果插值节点相同,在满足插值条件下用不同方法建立的插值公式是等价的【参考答案】: B26.在区间[xi-1,xi]上作线性插值,在图形上即为把两点用线段相连,n条线段组成折线,该折线对应的函数称为()A.牛顿插值函数B.分段线性插值函数C.三次样条插值函数D.拉格朗日插值函数【参考答案】: B27.A.AB.BC.CD.D【参考答案】: C28.牛顿插值多项式的优点是在计算时,高一级的插值多项式可利用前一次插值的结果A.正确B.错误【参考答案】: A29.f(x)=x^2+1,则f[1,2,3,4]=()A.4B.3C.1D.0【参考答案】: D30.A.AB.BC.CD.D 【参考答案】: C31.A.AB.BC.CD.D 【参考答案】: D32.有数表x=0,0.5,1,1.5,2,2.5;f(x)=-2,-1.75,-1,0.25,2,4.25.所确定的插值多项式次数是()A.二次B.三次C.四次D.五次【参考答案】: A33.若方阵A的谱半径小于1,则解方程组Ax=b的Jacobi迭代法收敛A.正确B.错误【参考答案】: B34.A.AB.BC.CD.D【参考答案】: B35.用二分法求方程x3+x-4=0在区间[1,2]内的根,要求精确到第三位小数,需对分()次A.8B.11C.9D.10【参考答案】: D36.若线性方程组Ax=b的系数矩阵A为严格对角占优矩阵,则解方程组的Jacobi迭代法和Gauss-Seidel迭代法()A.都收敛B.都发散C.Jacobi迭代法收敛,Gauss-Seidel迭代法发散 D.Jacobi迭代法发散,Gauss-Seidel迭代法收敛【参考答案】: A37.用简单迭代法求方程f(x)=0的实根,把方程f(x)=0 表示成x=j(x),则f(x)的根是()A.y=j(x)与x轴交点的横坐标B.y=j(x)与y=x交点的横坐标C.y=x与x 轴的交点的横坐标D.y=j(x)与y=x交点【参考答案】: B38.A.AB.BC.CD.D【参考答案】: C39.A.AB.BC.CD.D【参考答案】: B40.一阶均差f(x0,x1)=(f(x0)-f(x1))/(x0-x1)A.正确B.错误【参考答案】: A。
数值计算大作业题目一、非线性方程求根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以上表示收入差距悬殊。
第一章 绪论1.把下列各数按四舍五入规则舍入为有3位小数的近似数,并写出近似数的绝对误差和相对误差,指出近似数有几位有效数字: 93.1822 4.32250 15.9477 17.3675 2.按四舍五入原则,将下列数舍成五位有效数字:816.9567 6.000015 17.32250 1.235651 93.18213 0.015236233.设 **,671.3,6716.3x x x 则==有几位有效数字? 4.若1/4用0.25来表示,问有多少位有效数字? 5.若 1.1062,0.947a b ==是经过舍入后得到的近似值,问:,a b a b +⨯各有几位有效数字?6.设120.9863,0.0062y y ==是经过舍入后作为12,x x 的近似值, 求11y 和21y 的计算值与真值的相对误差限及12y y 和得到真值的相对误差限. 7.设0,x x >的相对误差为δ,求ln x 的绝对误差.8.正方形的边长约为100cm ,应该怎样测量,才能使其面积的误差不超过12cm . 9.设x 的相对误差为a %,求x n 的相对误差.10.计算球的体积,为了使相对误差限为1%,问度量半径R 时允许的相对误差限如何?11.5631.2*=x 是经四舍五入得到的近似值,则其相对误差≤*r e __________ 12. 设 0000073.0,1416.3,1415926.3**=-==x x x x 则称_________误差13.设⎰+=1061dx xx I nn ,设计一个计算10I 的算法,并说明你的算法的合理性。
14.设028Y =,按递推公式1n n Y Y -= (1,2,n = ), 计算到100Y27.982≈(5位有效数字),试问计算100Y 将有多大误差. 15.求方程25610x x -+=的两个根,使它至少具有4位有效数字27.982≈).16.当N 充分大时,怎样求121d 1N N x x ++⎰?17.序列{}n y 满足递推关系101n n y y =- (1,2,n = ),若0 1.41y =≈(三位有效数字),计算到10y 时误差有多大?这个计算过程稳定吗?18.计算61)f =1.41≈,利用下列算式计算,哪一个得到的结果最好?,3(3-,99- 19.()ln(f x x =,求(30)f 的值,若开平方用6位函数表,问求对数时误差有多大?若改用另一等价公式ln(ln(x x =- 计算,求对数时误差有多大?第二章 解线性方程组的直接方法1.用高斯消去法解方程组123234011921261x x x ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥-⎣⎦⎣⎦⎣⎦ 2.用LU 分解,将上题系数矩阵分解为L 和U 的乘积,L 是对角线元素为1的下三角矩阵,U 是上三角矩阵。
习题一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 )1cos2xy x-=;(4)(A)9y =-(B )y =解:当两个同(异)号相近数相减(加)时,相对误差可能很大,会严重丧失有效数字;当两个数相乘(除)时,大因子(小除数)可能使积(商)的绝对值误差增大许多。
课程设计课程名称:设计题目:学号:姓名:完成时间:题目一:非线性方程求根 一 摘要非线性方程的解析解通常很难给出,因此非线性方程的数值解就尤为重要。
本实验通过使用常用的求解方法二分法和Newton 法及改进的Newton 法处理几个题目,分析并总结不同方法处理问题的优缺点。
观察迭代次数,收敛速度及初值选取对迭代的影响。
用Newton 法计算下列方程(1) 310x x --= , 初值分别为01x =,00.45x =,00.65x =; (2) 32943892940x x x +-+= 其三个根分别为1,3,98-。
当选择初值02x =时给出结果并分析现象,当6510ε-=⨯,迭代停止。
解:1)采用MATLAB 进行计算;首先定义了Newton 法:function kk=newton(f,df,x0,tol,N)% Newton Method (牛顿法)% The first parameter f is a external function with respect to viable x.(第一个参数也就是本题所用的函数f )% The second parameter df is the first order diffential function of fx.(第二个参数也就是本体所用函数f 的导数方程df ) % x0 is initial iteration point(初值). % tol is the tolerance of the loop (精度).% N is the maximum number of iterations (循环上限). x=x0;f0=eval(f);df0=eval(df); n=0;disp(' [ n xn xn+1 fn+1 ]'); while n<=N x1=x0-f0/df0; x=x1; f1=eval(f); X=[n,x0,x1,f1]; disp(X);if abs(x0-x1)<tolfprintf('The procedure was successful.') kk=X; return else n=n+1; x0=x1;f0=f1;endendif n==N+1fprintf('the method failed after N iterations. '),kk=0;End我们把Newton法存为.m格式的文件;之后我们运行程序:clear;clc;syms xf=x^3-x-1;df=diff(f,x);x=newton(f,df,1,0.0001,50);x会得到一下结果[ n xn xn+1 fn+1 ]0 1.0000 1.5000 0.87501.0000 1.5000 1.0625 -0.86302.0000 1.0625 1.4940 0.8408到第50次迭代时候会出现该问题:47.0000 1.4898 1.0814 -0.816748.0000 1.0814 1.4898 0.816749.0000 1.4898 1.0814 -0.816750.0000 1.0814 1.4898 0.8167the method failed after N iterations.x =0;同样测试x0=0.45、0.65得不出结果,判断出初值离真值太远,所以我们采用牛顿下山法进行计算迭代:我们定义了其中的f函数和df函数,并且分别存为.m格式的文件,其代码如下:f:function y=f(x)y=x^3-x-1;df:function y=df(x)y=3*x^2-1;之后我们定义newton下山法同时也存为.m的程序:function [x,i]=downnewton(f,df,x0,tol)k=0;i=1;disp(' [ n xn xn+1 fn+1 ]'); while(k==0)fx=feval('f',x0);dfx=feval('df',x0);t=0;u=1;while(t==0)dx=-fx/dfx;x1=x0+u*dx;fx1=feval('f',x1);fx0=feval('f',x0);if(abs(fx1)>abs(fx0));u=u/2;elset=1;endendX=[i,x0,x1,fx1];disp(X);if(abs(fx1)<tol)k=1;elsex0=x1;i=i+1;endendx=x1;i=i;end之后带入x0=0.45;downnewton('f','df',0.45,10^(-6))[ n xn xn+1 fn+1 ]1.0000 0.4500 -0.4155 -0.65622.0000 -0.4155 -0.5857 -0.61523.0000 -0.5857 -0.5754 -0.61514.0000 -0.5754 -0.5782 -0.61515.0000 -0.5782 -0.5773 -0.61516.0000 -0.5773 -0.5774 -0.61517.0000 -0.5774 -0.5773 -0.61518.0000 -0.5773 -0.5774 -0.61519.0000 -0.5774 -0.5774 -0.615110.0000 -0.5774 -0.5774 -0.615111.0000 -0.5774 1.3131 -0.049012.0000 1.3131 1.3248 0.000513.0000 1.3248 1.3247 0.0000ans =1.3247带入x0=0.6;downnewton('f','df',0.6,10^(-6))[ n xn xn+1 fn+1 ]1.0000 0.6000 1.1406 -0.65662.0000 1.1406 1.3668 0.18663.0000 1.3668 1.3263 0.00674.0000 1.3263 1.3247 0.00005.0000 1.3247 1.3247 0.0000ans =1.3247带入x0=1;downnewton('f','df',1,10^(-6))[ n xn xn+1 fn+1 ]1.0000 1.0000 1.5000 0.87502.0000 1.5000 1.3478 0.10073.0000 1.3478 1.3252 0.00214.0000 1.3252 1.3247 0.0000ans =1.32472)同样采用Newton下山法:重新定义f、df:f:function y=f(x)y=x^3+94*x^2-389*x+294;df:function y=df(x)y=3*x^2+188*x-389;再带入初值x0=2;downnewton('f','df',2,5*10^(-6))[ n xn xn+1 fn+1 ]1 2 -98 0ans =-98得出x=-98;分析:先画出该函数的图像;x=(-100:.1:100);ezplot('x^3+94*x^2-389*x+294',[-100 100]) 得出该图像如图:-100-80-60-40-20020406080024681012141618x 105xx 3+94 x 2-389 x+294根据牛顿法的几何解释,在x0=2的点做切线,与y 相交,交点的横坐标值为x=-98则结束了该现象。
数值分析上机作业第 1 章1.1计算积分,n=9。
(要求计算结果具有6位有效数字)程序:n=1:19;I=zeros(1,19);I(19)=1/2*((exp(-1)/20)+(1/20));I(18)=1/2*((exp(-1)/19)+(1/19));for i=2:10I(19-i)=1/(20-i)*(1-I(20-i));endformat longdisp(I(1:19))结果截图及分析:在MATLAB中运行以上代码,得到结果如下图所示:当计算到数列的第10项时,所得的结果即为n=9时的准确积分值。
取6位有效数字可得.1.2分别将区间[-10.10]分为100,200,400等份,利用mesh或surf命令画出二元函数z=的三维图形。
程序:>> x = -10:0.1:10;y = -10:0.1:10;[X,Y] = meshgrid(x,y);Z = exp(-abs(X))+cos(X+Y)+1./(X.^2+Y.^2+1);subplot(2,2,1);mesh(X,Y,Z);title('步长0.1')>> x = -10:0.2:10;y = -10:0.2:10;[X,Y] = meshgrid(x,y);Z = exp(-abs(X))+cos(X+Y)+1./(X.^2+Y.^2+1);subplot(2,2,1);mesh(X,Y,Z);title('步长 0.2')>>x = -10:0.05:10;y = -10:0.05:10;[X,Y] = meshgrid(x,y);Z = exp(-abs(X))+cos(X+Y)+1./(X.^2+Y.^2+1);subplot(2,2,1);mesh(X,Y,Z);title('步长0.05')结果截图及分析:由图可知,步长越小时,绘得的图形越精确。
数值计算方法大作业目录第一章非线性方程求根 (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处的值。