北航数值分析计算实习报告一

合集下载

北航数值分析实验报告

北航数值分析实验报告

北航‎数值‎分析‎实验‎报告‎‎篇一‎:‎北航‎数值‎分析‎报告‎第一‎大题‎《‎数值‎分析‎》计‎算实‎习报‎告‎第一‎大题‎学‎号:‎D‎Y1‎30‎5‎姓名‎:‎指导‎老师‎:‎一、‎题目‎要求‎已‎知5‎01‎*5‎01‎阶的‎带状‎矩阵‎A,‎其特‎征值‎满足‎?1‎?‎2‎..‎.‎?5‎01‎。

试‎求:‎1‎、?‎1,‎?5‎01‎和?‎s的‎值;‎‎2、‎A的‎与数‎?k‎??‎1?‎k‎?5‎01‎??‎1‎40‎最‎接近‎的特‎征值‎?i‎k(‎k=‎1,‎2,‎..‎.,‎39‎);‎‎3、‎A的‎(谱‎范数‎)条‎件数‎c n‎d(‎A)‎2和‎行列‎式d‎e t‎A。

‎‎二、‎算法‎设计‎方案‎题‎目所‎给的‎矩阵‎阶数‎过大‎,必‎须经‎过去‎零压‎缩后‎进行‎存储‎和运‎算,‎本算‎法中‎压缩‎后的‎矩阵‎A1‎如下‎所示‎。

‎?0‎?0‎?A‎1?‎?a‎1‎??‎b?‎?c‎0‎b a‎2b‎c‎c b‎b c‎.‎..‎..‎..‎..‎..‎.‎c b‎b c‎c‎b a‎50‎0b‎0‎a ‎3.‎..‎a4‎99‎c‎?‎b?‎?a‎50‎1?‎?‎0?‎0?‎?‎由矩‎阵A‎的特‎征值‎满足‎的条‎件可‎知‎?1‎与?‎50‎1之‎间必‎有一‎个最‎大,‎则采‎用幂‎法求‎出的‎一‎个特‎征值‎必为‎其中‎的一‎个:‎当‎所求‎得的‎特征‎值为‎正数‎,则‎为?‎50‎1;‎否则‎为?‎1。

‎在求‎得?‎1与‎?‎50‎1其‎中的‎一个‎后,‎采用‎带位‎移的‎幂法‎则可‎求出‎它们‎中的‎另一‎个,‎且位‎移量‎即为‎先求‎出的‎特‎征值‎的值‎。

用‎反幂‎法求‎得的‎特征‎值必‎为?‎s。

‎由条‎件数‎的性‎质可‎得,‎c n‎d(‎A)‎2为‎模最‎大的‎特征‎值与‎模最‎小的‎特征‎值之‎比的‎模,‎因此‎,求‎出?‎1,‎?5‎01‎和?‎s的‎值后‎,则‎可以‎求得‎c n‎d(‎A)‎2。

北航数值分析计算实习1

北航数值分析计算实习1

《数值分析》计算实习题目110091013 劳云杰一、算法设计方案根据提示的算法,首先使用幂法求出按模最大的特征值λt1,再根据已求出的λt1用带原点平移的幂法求出另一个特征值λt2,比较两个λ的大小,根据已知条件,可以得出λ1和λ501.至于λs,由于是按模最小的特征值,使用反幂法求之,由于反幂法需要解线性方程组,故对矩阵进行Doolittle分解。

再通过带原点平移的反幂法求跟矩阵的与数最接近的特征值。

对非奇异的矩阵A,根据条件数定义,取λt1/λs的绝对值,两个特征值在之前步骤中均以求得。

由于对矩阵进行了Doolittle分解,所以矩阵的行列式det A可由分解得出的上三角阵U 的对角线上元素相乘求得。

为了使A的所有零元素都不存储,使用书本25页的压缩存储法对A进行存储,在计算时通过函数在数组C中检索A中元素即可。

由于A是501*501矩阵,C应取为5*501矩阵。

由于数据不大,为了方便起见,在程序中取502*502矩阵或者502向量,C也取为6*502矩阵。

程序编写参考《数值分析》颜庆津著和[C数值算法].(美国)W ILLIAM.H.P RESS.扫描版。

二、全部源程序#include <stdio.h>#include <math.h>#define XS 1.0e-12//精度水平void fz_a();//对矩阵A赋值double js(int,int);//在压缩矩阵中检索A的元素double mf(double);//幂法double fmf(double);//反幂法int lu(double);//Doolittle分解int jfc(double[],double[]);//解方程int max(int,int);int min(int,int);double (*u)[502]=new double[502][502];//上三角阵double (*l)[502]=new double[502][502];//单位下三角阵double a[6][502];//压缩存储矩阵int max(int x,int y)//比大小函数×2{ return (x>y?x:y);}int min(int x,int y)//精度关系,比较下标用{ return (x<y?x:y);}int main(){printf("请耐心等待,先看看中间过程吧~\n");int i,k;double ldt1,ldt2,ld1,ld501,lds,mu[40],det;double ld[40];fz_a();//对A赋值ldt1=mf(0);//幂法求模最大的特征值ldt2=mf(ldt1);//以第一次求得的特征值进行平移ld1=ldt1<ldt2?ldt1:ldt2;//大的就是λ501ld501=ldt1<ldt2?ldt2:ldt1;lu(0);lds=fmf(0);//反幂法求λsdet=1;//初始化行列式for(i=1;i<=501;i++)det=det*u[i][i];//用U的对角元素求行列式for(k=1;k<=39;k++){mu[k]=ld1+k*(ld501-ld1)/40;//与数lu(mu[k]);ld[k]=fmf(mu[k]);}printf("\n 列出结果\n");printf("λ1=%1.12e λ501=%1.12e\n",ld1,ld501);printf("λs=%1.12e \n",lds);printf("cond(A)=%1.12e \n",fabs(ldt1/lds));printf("detA=%1.12e \n",det);for(k=1;k<=39;k++)//列出跟与数最接近特征值{printf("λi%d=%1.12e\t",k,ld[k]);if(k%2==0)printf("\n");}//界面友好性delete []u;delete []l;getchar();return 0;}void fz_a()//对A赋值{int i;for(i=3;i<=501;i++)a[1][i]=a[5][502-i]=-0.064;//原A矩阵的cfor(i=2;i<=501;i++)a[2][i]=a[4][502-i]=0.16;//原A矩阵的bfor(i=1;i<=501;i++)a[3][i]=(1.64-0.024*i)*sin(0.2*i)-0.64*exp(0.1/i);//原对角线元素}double js(int i,int j)//对压缩矩阵检索A的元素{if(abs(i-j)<=2)return a[i-j+3][j];else return 0;}double mf(double offset)//幂法{int i,x1;double u[502],y[502];double beta=0,prebeta=-1000,yita=0;//用幂法的第一种迭代方法for(i=1;i<=501;i++) //用到了2-范数u[i]=1,y[i]=0;for(int k=1;k<=10000;k++)//对迭代次数进行限制{yita=0;for(i=1;i<=501;i++)yita=sqrt(yita*yita+u[i]*u[i]);for(i=1;i<=501;i++)y[i]=u[i]/yita;for(x1=1;x1<=501;x1++){u[x1]=0;for(int x2=1;x2<=501;x2++)u[x1]=u[x1]+((x1==x2)?(js(x1,x2)-offset):js(x1,x2))*y[x2];}prebeta=beta;beta=0;for(i=1;i<=501;i++)beta=beta+y[i]*u[i];if(fabs((prebeta-beta)/beta)<=XS){printf("offset=%f lb=%f err=%e k=%d\n",offset,(beta+offset),fabs((prebeta-beta)/beta),k);break;};}//满足误差条件后,迭代终止,并输出平移量,误差和迭代次数return(beta+offset);//加上平移量,方便比较}double fmf(double offset)//反幂法{ int i;double u[502],y[502];double beta=0,prebeta=0,yita=0;for(i=1;i<=501;i++)u[i]=1,y[i]=0; //相关量初始化for(int k=1;k<=10000;k++)//限制迭代次数{yita=0;for(i=1;i<=501;i++)yita=sqrt(yita*yita+u[i]*u[i]);for(i=1;i<=501;i++)y[i]=u[i]/yita;jfc(u,y);prebeta=beta;beta=0;for(i=1;i<=501;i++)beta=beta+y[i]*u[i];beta=1/beta;if(fabs((prebeta-beta)/beta)<=XS){printf("offset=%f lb=%f err=%ek=%d\n",offset,(beta+offset),fabs((prebeta-beta)/beta),k);break;};}//满足误差条件后,迭代终止,并输出平移量,误差和迭代次数return(beta+offset);}int lu(double offset)//Doolittle分解{int i,j,k,t;double sum;//中间量for(k=1;k<=501;k++)for(j=1;j<=501;j++){u[k][j]=l[k][j]=0;if(k==j)l[k][j]=1;}//对LU矩阵初始化for(k=1;k<=501;k++)//对式(2.12)的程序实现{for(j=k;j<=min(k+2,501);j++){sum=0;for(t=max(1,max(k-2,j-2));t<=(k-1);t++)sum=sum+l[k][t]*u[t][j];//j=k,k+1,……,nu[k][j]=((k==j)?(js(k,j)-offset):js(k,j))-sum;}if(k==501)continue;for(i=k+1;i<=min(k+2,501);i++)//i=k+1,……,n{sum=0;for(t=max(1,max(i-2,k-2));t<=(k-1);t++)sum=sum+l[i][t]*u[t][k];l[i][k]=(((i==k)?(js(i,k)-offset):js(i,k))-sum)/u[k][k];}}return 0;}int jfc(double x[],double b[])//解方程{int i,t;double y[502];double sum;y[1]=b[1];for(i=2;i<=501;i++){sum=0;for(t=max(1,i-2);t<=i-1;t++)sum=sum+l[i][t]*y[t];y[i]=b[i]-sum;}x[501]=y[501]/u[501][501];for(i=500;i>=1;i--){sum=0;for(t=i+1;t<=min(i+2,501);t++)sum=sum+u[i][t]*x[t];x[i]=(y[i]-sum)/u[i][i];}return 0;}三、结果λ1=-1.070011361502e+001λ501=9.724634098777e+000λs=-5.557910794230e-003cond(A)=1.925204273902e+003detA=2.772786141752e+118λi1=-1.018293403315e+001 λi2=-9.585707425068e+000 λi3=-9.172672423928e+000λi4=-8.652284007898e+000 λi5=-8.0934********e+000 λi6=-7.659405407692e+000λi7=-7.119684648691e+000 λi8=-6.611764339397e+000 λi9=-6.0661********e+000λi10=-5.585101052628e+000 λi11=-5.114083529812e+000 λi12=-4.578872176865e+000λi13=-4.096470926260e+000 λi14=-3.554211215751e+000 λi15=-3.0410********e+000 λi16=-2.533970311130e+000 λi17=-2.003230769563e+000 λi18=-1.503557611227e+000 λi19=-9.935586060075e -001 λi20=-4.870426738850e -001 λi21=2.231736249575e -002 λi22=5.324174742069e -001 λi23=1.052898962693e+000 λi24=1.589445881881e+000 λi25=2.060330460274e+000 λi26=2.558075597073e+000 λi27=3.080240509307e+000 λi28=3.613620867692e+000 λi29=4.0913********e+000 λi30=4.603035378279e+000 λi31=5.132924283898e+000 λi32=5.594906348083e+000 λi33=6.080933857027e+000 λi34=6.680354092112e+000 λi35=7.293877448127e+000 λi36=7.717111714236e+000 λi37=8.225220014050e+000 λi38=8.648666065193e+000 λi39=9.254200344575e+000四、讨论迭代初始向量的选取对计算结果的影响1.在反幂法中取迭代向量u[1]=1,u[i]=0,i=2,……,501,最后得出的结果中λs=2.668886923785e -002,cond(A)也随之改变成4.009204556274e+0022.在幂法中取迭代向量u[1]=1,u[i]=2,i=2,……,501,最后得出的结果不变。

数值分析计算方法实验报告

数值分析计算方法实验报告
break;
end;
end;
X=x;
disp('迭代结果:');
X
format short;
输出结果:
因为不收敛,故出现上述情况。
4.超松弛迭代法:
%SOR法求解实验1
%w=1.45
%方程组系数矩阵
clc;
A=[2,10,0,-3;-3,-4,-12,13;1,2,3,-4;4,14,9,-13]
b=[10,5,-2,7]'
b=[10,5,-2,7]'
[m,n]=size(A);
if m~=n
error('矩阵A的行数和列数必须相同');
return;
end
if m~=size(b)
error('b的大小必须和A的行数或A的列数相同');
return;
end
if rank(A)~=rank([A,b])
error('A矩阵的秩和增广矩阵的秩不相同,方程不存在唯一解');
3.实验环境及实验文件存档名
写出实验环境及实验文件存档名
4.实验结果及分析
输出计算结果,结果分析和小结等。
解:1.高斯列主元消去法:
%用高斯列主元消去法解实验1
%高斯列主元消元法求解线性方程组Ax=b
%A为输入矩阵系数,b为方程组右端系数
%方程组的解保存在x变量中
format long;
A=[2,10,0,-3;-3,-4,-12,13;1,2,3,-4;4,14,9,-13]
return;
end
c=n+1;
A(:,c)=b;
for k=1:n-1

计算实习报告4篇【实用模板】

计算实习报告4篇【实用模板】

计算实习报告4篇计算实习报告篇1(2110字)1、实习目的通过理论联系实际,巩固所学的知识,提高处理实际问题的能力,为顺利毕业进行做好充分的准备,并为自己能顺利与社会环境接轨做准备。

通过这次实习,使我们进一步理解和领会所学的基本理论,了解计算机技术和信息管理技术的发展及应用,较为系统地掌握计算机应用技能和信息管理技能,把所学知识与解决实际问题相联系,能够利用计算机处理工作中的各种信息,培养我们发现问题、分析问题和解决问题的能力,从而提高我们从事实际工作的能力。

2、实习意义生产实习是一个极为重要的实践性教学环节。

通过实习,使学生在社会实践中接触与本专业相关的实际工作,增强感性认识,培养和锻炼学生综合运用所学的基础理论、基本技能和专业知识,去独立分析和解决实际问题的能力,把理论和实践结合起来,提高实践动手能力,为学生毕业后走上工作岗位打下一定的基础;同时可以检验教学效果,为进一步提高教育教学质量,培养合格人才积累经验。

计算机是一门对实践要求较高的学科,通过专业实习,使学生能熟悉有关计算机专业的各个领域,使学生毕业后能胜任与本专业相关的工作。

大学继续教育5年中学习了很多,经历了很多,得到的是学习能力、处事能力和一些专业知识。

可面对社会,我们经验太少,思想单纯!毕业实习,给了我们一个了解社会,增加经验,熟悉工作单位的机会。

锻炼自己的动手能力,将学习的理论知识运用于实践当中,反过来还能检验书本上理论的正确性,有利于融会贯通。

同时,也能开拓视野,完善自己的知识结构,达到锻炼能力的目的。

一切都是为了让实践者对本专业知识形成一个客观,理性的认识,从而不与社会现实相脱节。

此外通过理论联系实际,巩固所学的知识,提高处理实际问题的能力,了解设计专题的主要内容,为毕业设计的顺利进行做好充分的准备,并为自己能顺利与社会环境接轨做准备。

3、实习单位调研情况我的实习单位中广有线信息网络有限公司枣庄分公司(以下简称中广有线枣庄分公司)是枣庄市广播电视局国有资产出资参股成立的有限责任企业,坐落在美丽的光明广场西侧。

数值分析实验报告心得(3篇)

数值分析实验报告心得(3篇)

第1篇在数值分析这门课程的学习过程中,我深刻体会到了理论知识与实践操作相结合的重要性。

通过一系列的实验,我对数值分析的基本概念、方法和应用有了更加深入的理解。

以下是我对数值分析实验的心得体会。

一、实验目的与意义1. 巩固数值分析理论知识:通过实验,将课堂上学到的理论知识应用到实际问题中,加深对数值分析概念和方法的理解。

2. 培养实际操作能力:实验过程中,我学会了使用Matlab等软件进行数值计算,提高了编程能力。

3. 增强解决实际问题的能力:实验项目涉及多个领域,通过解决实际问题,提高了我的问题分析和解决能力。

4. 培养团队协作精神:实验过程中,我与同学们分工合作,共同完成任务,培养了团队协作精神。

二、实验内容及方法1. 实验一:拉格朗日插值法与牛顿插值法(1)实验目的:掌握拉格朗日插值法和牛顿插值法的原理,能够运用这两种方法进行函数逼近。

(2)实验方法:首先,我们选择一组数据点,然后利用拉格朗日插值法和牛顿插值法构造插值多项式。

最后,我们将插值多项式与原始函数进行比较,分析误差。

2. 实验二:方程求根(1)实验目的:掌握二分法、Newton法、不动点迭代法、弦截法等方程求根方法,能够运用这些方法求解非线性方程的根。

(2)实验方法:首先,我们选择一个非线性方程,然后运用二分法、Newton法、不动点迭代法、弦截法等方法求解方程的根。

最后,比较不同方法的收敛速度和精度。

3. 实验三:线性方程组求解(1)实验目的:掌握高斯消元法、矩阵分解法等线性方程组求解方法,能够运用这些方法求解线性方程组。

(2)实验方法:首先,我们构造一个线性方程组,然后运用高斯消元法、矩阵分解法等方法求解方程组。

最后,比较不同方法的计算量和精度。

4. 实验四:多元统计分析(1)实验目的:掌握多元统计分析的基本方法,能够运用这些方法对数据进行分析。

(2)实验方法:首先,我们收集一组多元数据,然后运用主成分分析、因子分析等方法对数据进行降维。

数值分析实验报告5篇

数值分析实验报告5篇
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -11 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -12 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -13 -14
1.69376699767424 0.92310666706964 0.08471614569741 0.40804026409411
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
讨论:
利用这种方法进行这类实验,可以很精确的扰动敏感性的一般规律。即 当对扰动项的系数越来越小时,对其多项式扰动的结果也就越来越小, 即扰动敏感性与扰动项的系数成正比,扰动项的系数越大,对其根的扰 动敏感性就越明显,当扰动的系数一定时,扰动敏感性与扰动的项的幂 数成正比,扰动的项的幂数越高,对其根的扰动敏感性就越明显。
解线性方程组的直接方法
实验 (主元的选取与算法的稳定性) 问题提出:Gauss消去法是我们在线性代数中已经熟悉的。但由于计算 机的数值运算是在一个有限的浮点数集合上进行的,如何才能确保 Gauss消去法作为数值算法的稳定性呢?Gauss消去法从理论算法到数值 算法,其关键是主元的选择。主元的选择从数学理论上看起来平凡,它 却是数值分析中十分典型的问题。 实验内容:考虑线性方程组 编制一个能自动选取主元,又能手动选取主元的求解线性方程组的 Gauss消去过程。 实验要求: (1)取矩阵,则方程有解。取n=10计算矩阵的条件数。让程序自动选 取主元,结果如何? (2)现选择程序中手动选取主元的功能。每步消去过程总选取按模最 小或按模尽可能小的元素作为主元,观察并记录计算结果。若每步消去 过程总选取按模最大的元素作为主元,结果又如何?分析实验的结果。 (3)取矩阵阶数n=20或者更大,重复上述实验过程,观察记录并分析 不同的问题及消去过程中选择不同的主元时计算结果的差异,说明主元

数值分析实验报告总结

数值分析实验报告总结

一、实验背景数值分析是研究数值计算方法及其理论的学科,是计算机科学、数学、物理学等领域的重要基础。

为了提高自身对数值分析理论和方法的理解,我们进行了数值分析实验,通过实验加深对理论知识的掌握,提高实际操作能力。

二、实验目的1. 理解数值分析的基本理论和方法;2. 掌握数值分析实验的基本步骤和技巧;3. 培养实验设计和数据分析能力;4. 提高编程和计算能力。

三、实验内容本次实验主要分为以下几个部分:1. 线性方程组求解实验:通过高斯消元法、LU分解法等求解线性方程组,并分析算法的稳定性和误差;2. 矩阵特征值问题计算实验:利用幂法、逆幂法等计算矩阵的特征值和特征向量,分析算法的收敛性和精度;3. 非线性方程求根实验:运用二分法、牛顿法、不动点迭代法等求解非线性方程的根,比较不同算法的优缺点;4. 函数插值实验:运用拉格朗日插值、牛顿插值等方法对给定的函数进行插值,分析插值误差;5. 常微分方程初值问题数值解法实验:运用欧拉法、改进的欧拉法、龙格-库塔法等求解常微分方程初值问题,比较不同算法的稳定性和精度。

四、实验过程1. 线性方程组求解实验:首先,编写程序实现高斯消元法、LU分解法等算法;然后,对给定的线性方程组进行求解,记录计算结果;最后,分析算法的稳定性和误差。

2. 矩阵特征值问题计算实验:编写程序实现幂法、逆幂法等算法;然后,对给定的矩阵进行特征值和特征向量的计算,记录计算结果;最后,分析算法的收敛性和精度。

3. 非线性方程求根实验:编写程序实现二分法、牛顿法、不动点迭代法等算法;然后,对给定的非线性方程进行求根,记录计算结果;最后,比较不同算法的优缺点。

4. 函数插值实验:编写程序实现拉格朗日插值、牛顿插值等方法;然后,对给定的函数进行插值,记录计算结果;最后,分析插值误差。

5. 常微分方程初值问题数值解法实验:编写程序实现欧拉法、改进的欧拉法、龙格-库塔法等算法;然后,对给定的常微分方程初值问题进行求解,记录计算结果;最后,比较不同算法的稳定性和精度。

数值分析 实验报告

数值分析 实验报告

数值分析实验报告1. 引言数值分析是一门研究如何利用计算机进行数值计算的学科。

它涵盖了数值计算方法、数值逼近、插值和拟合、数值微积分等内容。

本实验报告旨在介绍数值分析的基本概念,并通过实验验证其中一些常用的数值计算方法的准确性和可行性。

2. 实验目的本实验的目的是通过对一些简单数学问题进行数值计算,验证数值计算方法的正确性,并分析计算误差。

具体实验目标包括: - 了解数值计算方法的基本原理和应用场景; - 掌握常用的数值计算方法,如二分法、牛顿法等; - 验证数值计算方法的准确性和可靠性。

3. 实验设计3.1 实验问题选择了以下两个数学问题作为实验对象: 1. 求解方程f(x) = 0的根; 2. 求解函数f(x)在给定区间上的最小值。

3.2 实验步骤3.2.1 方程求根1.确定待求解的方程f(x) = 0;2.选择合适的数值计算方法,比如二分法、牛顿法等;3.编写相应的计算程序,并根据初始条件设置迭代终止条件;4.运行程序,得到方程的根,并计算误差。

3.2.2 函数最小值1.确定待求解的函数f(x)和给定的区间;2.选择合适的数值计算方法,比如黄金分割法、斐波那契法等;3.编写相应的计算程序,并根据初始条件设置迭代终止条件;4.运行程序,得到函数的最小值,并计算误差。

4. 实验结果与分析4.1 方程求根我们选择了二分法和牛顿法来求解方程f(x) = 0的根,并得到了如下结果: - 二分法得到的根为 x = 2.345,误差为 0.001; - 牛顿法得到的根为 x = 2.345,误差为 0.0001。

通过计算结果可以看出,二分法和牛顿法都能较准确地求得方程的根,并且牛顿法的收敛速度更快。

4.2 函数最小值我们选择了黄金分割法和斐波那契法来求解函数f(x)在给定区间上的最小值,并得到了如下结果: - 黄金分割法得到的最小值为 x = 3.142,误差为 0.001; - 斐波那契法得到的最小值为 x = 3.142,误差为 0.0001。

北航数值分析计算实习第一题编程

北航数值分析计算实习第一题编程

i − t + s +1,t t − k + s +1, k t = max(1,i − r ,k − s )
∑c
c
) / cs +1, k
[i = k + 1, k + 2,⋯ , min( k + r , n); k < n]
(2) 求解 Ly = b,Ux = y (数组 b 先是存放原方程右端向量,后来存放中间向量 y)
0 b a2
b c
c b a3 b c
⋯ ⋯ ⋯ ⋯ ⋯
c b a499 b c
c b a500 b 0
c ⎤ b ⎥ ⎥ a501 ⎥ ⎥ 0 ⎥ 0 ⎥ ⎦
在数组 C 中检索矩阵 A 的带内元素 aij 的方法是: A 的带内元素 aij =C 中的元素 ci − j + s +1, j
2
数值分析计算实习题目一
i −1
bi := bi −
பைடு நூலகம்
i − t + s +1,t t t = max(1,i − r )
∑c
b
(i = 2,3,⋯ , n)
xn := bn / cs +1, n
min( i + s )
xi := (bi −
t = i +1
∑c
i −t + s +1,t t
x ) / cs +1,i
(i = n − 1, n − 2,⋯ ,1)
3、Doolittle 分解求解 n 元带状线性方程组(doolittle()函数)
按照上述对带状矩阵 A 的存储方法和元素 aij 的检索方法,并且把三角分解的结果 ukj 和 lik 分 别存放在 akj 和 aik 原先的存储单元内,那么用 Doolittle 分解法求解 n 元带状线性方程组的算法 可重新表述如下(其中“:=”表示赋值) : (1) 作分解 A = LU 。 对于 k=1,2, ……,n 执行

数值分析实验报告1

数值分析实验报告1
end
p
得到m=(00)T
即M0=0 ;M1=;M2=;M3=;M4=0
则根据三次样条函数定义,可得:
S(x)=
接着,在Command Window里输入画图的程序代码,
下面是画牛顿插值以及三次样条插值图形的程序:
x=[ ];
y=[ ];
plot(x,y)
hold on
for i=1:1:5
y(i)= *(x(i)*(x(i)*(x(i)*(x(i)*(x(i)*(x(i)*(x(i)
Pn=f(x0)+f[x0,x1](x-x0)+ f[x0,x1,x2](x-x0) (x-x1)+···+ f[x0,x1,···xn](x-x0) ···(x-xn-1)
我们要知道牛顿插值多项式的系数,即均差表中得部分均差。
在MATLAB的Editor中输入程序代码,计算牛顿插值中多项式系数的程序如下:
【实验原理】
《数值分析》第二章“插值法”的相关内容,包括:牛顿多项式插值,三次样条插值,拉格朗日插值的相应算法和相关性质。
【实验环境】(使用的软硬件)
软件:
MATLAB 2012a
硬件:
电脑型号:联想 Lenovo 昭阳E46A笔记本电脑
操作系统:Windows 8 专业版
处理器:Intel(R)Core(TM)i3 CPU M 350 @
实验内容:
【实验方案设计】
第一步,将书上关于三种插值方法的内容转化成程序语言,用MATLAB实现;第二步,分别用牛顿多项式插值,三次样条插值,拉格朗日插值求解不同的问题。
【实验过程】(实验步骤、记录、数据、分析)
实验的主要步骤是:首先分析问题,根据分析设计MATLAB程序,利用程序算出问题答案,分析所得答案结果,再得出最后结论。

数值分析上机实习报告

数值分析上机实习报告

数值分析上机实习报告目录1.问题一 (1)问题一重述 (1)秦九韶算法简介 (1)问题一算法实现 (1)问题一求解 (1)2.问题二 (2)问题二重述 (2)逐次超松弛迭代法(SOR法)简介 (2)问题二算法实现 (3)问题二求解 (3)3.问题三 (4)问题三重述 (4)最小二乘拟合多项式与拉格朗日插值多项式简介 (4)3.2.1最小二乘拟合多项式简介 (4)3.2.2拉格朗日插值简介 (5)问题三算法实现 (5)3.3.1多项式拟合算法 (5)3.3.2拉格朗日插值算法 (6)问题三求解 (6)3.4.1最小二乘多项式拟合结果 (6)3.4.2拉格朗日插值结果 (8)问题三评判 (9)3.5.1问题三评判方式 (9)3.5.2问题三评判结果 (9)4.总结与体会 (10)5.附录 (11)1. 问题一问题一重述利用秦九韶算法简化求多项式1110n n n n x a x a y x a a --=++++的值的运算式,并写程序计算多项式42352x y x x =--+在1x =-点处的值。

秦九韶算法简介121210...n n n n y a x a x a x a x a --=+++++化为以下形式:1210(...(())...)n n n y a x a x a x a x a --=+++++求多项式值时先计算内层括号内的一次多项式的值,然后由内向外逐层计算一次多项式的值,即:11n n v a x a -=+212n v v x a -=+ …1k k n k v v x a +-=+…10n n v v x a -=+ 问题一算法实现Step1:输入多项式的降次排列的系数矩阵,某次缺失的系数用零补充之;Step2:计算表达式1v ,按递推1k k n k v v x a +-=+公式,一直计算到表达式n v ,表达式n v 即为所求秦九韶表达式;Step3:输入x 的值;Step4:计算1v ,按递推1k k n k v v x a +-=+公式,一直计算到n v 的值,n v 的值即为x 处多项式的值。

北航数值分析-实习作业1(C语言详细注释)

北航数值分析-实习作业1(C语言详细注释)

《数值分析》计算实习作业《一》北航第一题 设有501501⨯的矩阵⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=501500499321a bc b a b cc b a b ccb a bc c b a b c b a A其中.064.0,16.0);501,2,1(64.0)2.0sin()024.064.1(1.0-===--=c b i e i i a i i 矩阵的特征值)501,,2,1( =i i λ满足||min ||,501150121i i s λλλλλ≤≤=<<<试求1. 5011,λλ和s λ的值2. 的与数4015011λλκλμ-+=k 最接近的特征值)39,,2,1( =K κλi3. 的(谱范数)条件数2)A (cond 和行列式A det 要求1. 算法的设计方案(A 的所有零元素都不能存储)2. 全部源程序(详细注释)。

变量为double ,精度-1210=ε,输出为e 型12位有效数字3. 特征值s 5011,,λλλ和)39,,2,1( =K κλi 以及A cond det ,)A (2的值4. 讨论迭代初始向量的选取对计算结果的影响,并说明原因解答:1. 算法设计对于s λ满足||min ||5011i i s λλ≤≤=,所以s λ是按模最小的特征值,直接运用反幂法可求得。

对于5011,λλ,一个是最大的特征值,一个是最小的特征值,不能确定两者的绝对值是否相等,因此必须首先假设||||5011λλ≠,然后运用幂法,看能否求得一个特征值,如果可以求得一个,证明A 是收敛的,求得的结果是正确的,然后对A 进行带原点平移的幂法,偏移量是前面求得的特征值,可以求得另一个特征值,最后比较这两个特征值,较大的特征值是501λ,较小的特征值就是1λ。

如果在假设的前提下,无法运用幂法求得按模最大的特征值,即此时A 不收敛,则需要将A 进行带原点平移的幂法,平移量可以选取1,再重复上述步骤即可求得两个特征值。

2021年计算实习报告三篇

2021年计算实习报告三篇

2021年计算实习报告三篇计算实习报告篇1一、实习目的通过理论联系实际,巩固所学的知识,提高处理实际问题的能力,为顺利毕业进行做好充分的准备,并为自己能顺利与社会环境接轨做准备。

通过这次实习,使我们进一步理解和领会所学的基本理论,了解计算机技术和信息管理技术的发展及应用,较为系统地掌握计算机应用技能和信息管理技能,把所学知识与解决实际问题相联系,能够利用计算机处理工作中的各种信息,培养我们发现问题、分析问题和解决问题的能力,从而提高我们从事实际工作的能力。

通过理论联系实际,巩固所学的知识,提高处理实际问题的能力,了解设计专题的主要内容,使学生能够了解社会、学校的需要,在实习单位领导的帮助,对自己今后所从事的事业有一个实习了解的过程。

为毕业设计的顺利进行做好充分的准备,并为自己能顺利与社会环境接轨做准备。

二、实习意义通过实习,使学生在社会实践中接触与本专业相关的实际工作,增强感性认识,培养和锻炼学生综合运用所学的基础理论、基本技能和专业知识,去独立分析和解决实际问题的能力,把理论和实践结合起来,提高实践动手能力,为学生毕业后走上工作岗位打下一定的基础;同时可以检验教学效果,为进一步提高教育教学质量,培养合格人才积累经验。

计算机是一门对实践要求较高的学科,通过专业实习,使学生能熟悉有关计算机专业的各个领域,使学生毕业后能胜任与本专业相关的工作。

大学四年学习了很多,经历了很多,得到的是学习能力、处事能力和一些专业知识。

可面对社会,我们经验太少,思想单纯!毕业实习,给了我们一个了解社会,增加经验,熟悉工作单位的机会。

锻炼自己的动手能力,将学习的理论知识运用于实践当中,反过来还能检验书本上理论的正确性,有利于融会贯通。

同时,也能开拓视野,完善自己的知识结构,达到锻炼能力的目的。

一切都是为了让实践者对本专业知识形成一个客观,理性的认识,从而不与社会现实相脱节。

此外通过理论联系实际,巩固所学的知识,提高处理实际问题的能力,了解设计专题的主要内容,为毕业设计的顺利进行做好充分的准备,并为自己能顺利与社会环境接轨做准备。

北航数值分析计算实习报告一

北航数值分析计算实习报告一

北京航空航天大学《数值分析》计算实习报告第一大题学 院:自动化科学与电气工程学院 专 业: 控制科学与工程 学 生 姓 名: 学 号: 教 师: 电 话: 完 成 日 期: 2015年11月6日北京航空航天大学Beijing University of Aeronautics and Astronautics实习题目:第一题 设有501501⨯的实对称矩阵A ,其中,064.0,16.0),501,,2,1(64.0)2.0sin()024.064.1(1.0-==⋅⋅⋅=--=c b i e i i a ii 。

矩阵A 的特征值为)501,,2,1(⋅⋅⋅=i i λ,并且有1.求1λ,501λ和s λ的值。

2.求A 的与数4015011λλλμ-+=kk 最接近的特征值)39,,2,1(⋅⋅⋅=k k i λ。

3.求A 的(谱范数)条件数2)A (cond 和行列式detA 。

说明:1.在所用的算法中,凡是要给出精度水平ε的,都取12-10=ε。

2.选择算法时,应使矩阵A 的所有零元素都不储存。

3.打印以下内容: (1)全部源程序;(2)特征值),,39,...,2,1(,s 5011=k k i λλλλ以及A det ,)A (cond 2的值。

4.采用e 型输出实型数,并且至少显示12位有效数字。

一、算法设计方案1、求1λ,501λ和s λ的值。

由于||min ||,501150121i i s λλλλλ≤≤=≤⋅⋅⋅≤≤,可知绝对值最大特征值必为1λ和501λ其中之一,故可用幂法求出绝对值最大的特征值λ,如果λ=0,则1λ=λ,否则501λ=λ。

将矩阵A 进行一下平移:I -A A'λ= (1)对'A 用幂法求出其绝对值最大的特征值'λ,则A 的另一端点特征值1λ或501λ为'λ+λ。

s λ为按模最小特征值,||min ||5011i i s λλ≤≤=,可对A 使用反幂法求得。

数值分析实习报告

数值分析实习报告

数值分析实习报告姓名:***学号:***班级:***1.题目取h=0.1,利用Euler公式求解dy/dx=y-2*x/y (0<=x<=1)y(0)=12.思路利用左矩形公式得到的公式yn+1=yn+hf(xn,yn)进行迭代,在进行迭代的过程使用for循环,让n从1取到10,每迭代一次输出一个x与y,同时x加0.1再进行下一次迭代,一直到循环结束。

3.程序clear; y=1, x=0, %初始化for n=1:10y=1.1*y-0.2*x/y, x=x+0.1,end4.运行结果y =1 x =0y =1.1000 x =0.1000y =1.1918 x =0.2000y =1.2774 x =0.3000y =1.3582 x =0.4000y =1.4351 x = 0.5000y =1.5090 x = 0.6000y =1.5803 x = 0.7000y =1.6498 x = 0.8000y =1.7178 x =0.9000y =1.7848 x =1.00001.题目取h=0.1,利用Euler公式求解dy/dx=y-2*x/y (0<=x<=1)y(0)=12.思路利用右矩形公式得到的公式yn+1=yn+hf(xn+1,yn+1)进行迭代,在进行迭代的过程使用for 循环,让n从1取到10,每迭代一次输出一个x与y,同时x加0.1再进行下一次迭代,一直到循环结束。

3.程序clear; y0=1, x0=0, %初始化for n=1:10yp=y0+0.1*(y0-2*x0/y0);x=x0+0.1;x0=xyp=y0+0.1*(yp-2*x/yp)y0=ypend4.运行结果y =1 x =0y =1.0918 x =0.1000y =1.1763 x =0.2000y =1.2546 x =0.3000y =1.3278 x =0.4000y =1.3964 x =0.5000y =1.4609 x =0.6000y =1.5216 x =0.7000y =1.5786 x =0.8000y =1.6321 x =0.9000y =1.6819 x =1.00001.题目取h=0.1,利用梯形法求解dy/dx=y-2*x/y (0<=x<=1)y(0)=12.思路利用右矩形公式得到的公式yn+1=yn+hf(xn+1,yn+1)进行迭代,在进行迭代的过程使用for 循环,让n从1取到10,每迭代一次输出一个x与y,同时x加0.1再进行下一次迭代,一直到循环结束。

北航数值分析报告第一次大作业(幂法反幂法)

北航数值分析报告第一次大作业(幂法反幂法)

一、问题分析与算法描述1. 问题的提出:〔1〕用幂法、反幂法求矩阵的按摸最大和最小特征值,并求出相应的特征向量。

其中要求:迭代精度达到。

〔2〕用带双步位移的QR法求上述的全部特征值,并求出每一个实特征值相应的特征向量。

2. 算法的描述:(1) 幂法幂法主要用于计算矩阵的按摸为最大的特征值和相应的特征向量。

其迭代格式为:终止迭代的控制选用。

幂法的使用条件为实矩阵A具有n个线性无关的特征向量,其相应的特征值满足不等式或幂法收敛速度与比值或有关,比值越小,收敛速度越快。

(2) 反幂法反幂法用于计算实矩阵A按摸最小的特征值,其迭代格式为:每迭代一次都要求解一次线性方程组。

当k足够大时,,可近似的作为矩阵A的属于的特征向量。

比值越小,收敛的越快。

反幂法要求矩阵A非奇异。

(3) 带双步位移的QR分解法QR方法适用于计算一般实矩阵的全部特征值,尤其适用于计算中小型实矩阵的全部特征值。

本算例中采用带双步位移的QR方法,可加速收敛,其迭代格式为:二、计算结果与分析1. 计算结果:(1) 幂法:初始条件:最大迭代次数L=1000;向量计算结果:第1次迭代结果:最大特征值:0.00000e+000第2次迭代结果:最大特征值:2.48910e+000 相对误差:1.00000e+000 第3次迭代结果:最大特征值:1.67719e+000 相对误差:第4次迭代结果:最大特征值:-2.10960e+000 相对误差:1.79503e+000 第5次迭代结果:最大特征值:-6.13203e-001 相对误差:2.44030e+000 ……第794次迭代结果:最大特征值:-1.97638e+000 相对误差:最大特征值:-1.97638e+000 相对误差:********************最终迭代结果***************特征值:-1.97638e+000 相对误差:迭代次数:795(2) 反幂法:初始条件:最大迭代次数L=1000;向量运行结果:第1次迭代结果:最大特征值:1.07542e+000第2次迭代结果:最大特征值:-3.66550e+000 相对误差:1.29339e+000 第3次迭代结果:最大特征值:1.22709e+001 相对误差:1.29871e+000 第4次迭代结果:最大特征值:-1.03421e+000 相对误差:1.28650e+001 第5次迭代结果:最大特征值:相对误差:……第995次迭代结果:最大特征值:相对误差:第996次迭代结果:最大特征值:相对误差:最大特征值:相对误差:第998次迭代结果:最大特征值:相对误差:第999次迭代结果:最大特征值:相对误差:第1000次迭代结果:最大特征值:相对误差:******************************超过最大设定迭代次数,迭代失败!(3) 带双步位移的QR法:初始条件:最大迭代次数L=1000;向量运行结果:全部特征值:特征向量〔经谱X数归一化〕:实特征值对应特征向量:-0.062705 -0.022368 0.304372 0.064466 0.521833 -0.157024 0.136942 -0.218108 0.250264 -0.043064 -0.228688 -0.184632 -0.072871 0.124721 0.029070 0.102566 -0.136358 0.167727 0.085747 0.546165 实特征值对应特征向量:-0.018001 0.019652 0.273447 0.070528 0.274896 -0.144015 0.048385 0.376439 -0.583051 -0.054008 -0.168682 -0.113430 -0.034709 0.009204 0.472291 0.125664 -0.190617 0.113145 0.046278 0.059871 实特征值对应特征向量:0.106861 0.087709 -0.024967 -0.020897 0.064302 0.034047 0.535143 0.046383 0.028832 0.003479-0.097276 -0.383801 0.089445 -0.039560 -0.036928 -0.021330 0.014811 0.705836 -0.108904 0.082022 实特征值对应特征向量:-0.055201 0.003399 0.242191 0.102847 0.372470 -0.372826 0.113953 0.240659 -0.310401 -0.076590 -0.244632 -0.192549 -0.077259 0.263328 0.201662 0.154166 -0.407814 0.186782 0.094649 0.173302 实特征值对应特征向量:0.427828 -0.546801 0.007822 -0.382580 0.025199 0.012788 0.033241 0.005389 -0.004065 0.043524 -0.032112 -0.044233 0.135395 -0.006564 0.001214 0.020165 0.011678 0.050001 -0.585765 0.013115 实特征值对应特征向量:0.236032 -0.139250 -0.008143 0.638527 -0.009049 -0.002911 -0.001307 0.003054 0.006515 -0.030134 0.012712 0.011368 -0.018792 -0.001753 -0.005749 -0.014290 -0.005292 -0.014591 0.717590 0.001369 实特征值对应特征向量:-0.227404 -0.048154 0.022615 0.297305 0.070372 0.039927 0.078503 0.015822 -0.012182 0.605334 -0.083616 -0.106270 -0.573963 -0.019907 0.003839 0.051362 0.036567 0.115613 0.332707 0.036954 实特征值对应特征向量:-0.027768 -0.051081 -0.159642 -0.054573 -0.084441 0.118378 0.029553 0.211088 0.203867 0.0486272. 结果分析以上三种方法中,幂法计算共进展了795次迭代才达到收敛,计算量较大,收敛性不好;反幂法计算结果未能收敛,通过进一步分析发现,这是因为反幂法迭代程序未考虑按模最小特征值为复数的情况,造成迭代失败。

北航硕士研究生数值分析大作业一

北航硕士研究生数值分析大作业一

数值分析—计算实习作业一学院:17系专业:精密仪器及机械姓名:张大军学号:DY14171142014-11-11数值分析计算实现第一题报告一、算法方案算法方案如图1所示。

(此算法设计实现完全由本人独立完成)图1算法方案流程图二、全部源程序全部源程序如下所示#include <iostream.h>#include <iomanip.h>#include <math.h>int main(){double a[501];double vv[5][501];double d=0;double r[3];double uu;int i,k;double mifayunsuan(double *a,double weiyi);double fanmifayunsuan(double *a,double weiyi);void yasuo(double *A,double (*C)[501]);void LUfenjie(double (*C)[501]);//赋值语句for(i=1;i<=501;i++){a[i-1]=(1.64-0.024*i)*sin(0.2*i)-0.64*exp(0.1/i);}//程序一:使用幂方法求绝对值最大的特征值r[0]=mifayunsuan(a,d);//程序二:使用幂方法求求平移λ[0]后绝对值最大的λ,得到原矩阵中与最大特征值相距最远的特征值d=r[0];r[1]=mifayunsuan(a,d);//比较λ与λ-λ[0]的大小,由已知得if(r[0]>r[1]){d=r[0];r[0]=r[1];r[1]=d;}//程序三:使用反幂法求λr[2]=fanmifayunsuan(a,0);cout<<setiosflags(ios::right);cout<<"λ["<<1<<"]="<<setiosflags(ios::scientific)<<setprecision(12)<<r[0]<<endl;cout<<"λ["<<501<<"]="<<setiosflags(ios::scientific)<<setprecision(12)<<r[1]<<endl;cout<<"λ[s]="<<setiosflags(ios::scientific)<<setprecision(12)<<r[2]<<endl;//程序四:求A的与数u最接近的特征值for(k=1;k<40;k++){uu=r[0]+k*(r[1]-r[0])/40;cout<<"最接近u["<<k<<"]"<<"的特征值为"<<setiosflags(ios::scientific)<<setprecision(12)<<fanmifayunsuan(a,uu)<<endl;}//程序五:谱范数的条件数是绝对值最大的特征值除以绝对值最小的特征值的绝对值cout<<"cond(A)2="<<fabs(r[0]/r[2])<<endl;//程序六:A的行列式的值就是A分解成LU之U的对角线的乘积yasuo(a,vv);LUfenjie(vv);uu=1;for(i=0;i<501;i++){uu=uu*vv[2][i];}cout<<"Det(A)="<<uu<<endl;return 1;}double mifayunsuan(double *a,double weiyi){int i,k;double b=0.16;double c=-0.064;double ee,w,v1,v2,mm,sum;double u[501];double y[505]={0};for(i=0;i<501;i++)u[i]=1;//给u赋初值if (weiyi!=0){for (i=0;i<501;i++)a[i]-=weiyi;}ee=1;k=0;//使得初始计算时进入循环语句while(ee>1e-12){mm=0;for(i=0;i<501;i++){mm=mm+u[i]*u[i];}w=sqrt(mm);for(i=0;i<501;i++){y[i+2]=u[i]/w;//注意此处编程与书上不同,之后会解释它的巧妙之处1 }for(i=0;i<501;i++){u[i]=c*y[i]+b*y[i+1]+a[i]*y[i+2]+b*y[i+3]+c*y[i+4];//1显然巧妙之处凸显出来}sum=0;for(i=0;i<501;i++){sum+=y[i+2]*u[i];}v1=v2;v2=sum;//去除特殊情况,减少漏洞if(k==0){k++;}else{ee=fabs(v2-v1)/fabs(v2);}}if (weiyi!=0){for (i=0;i<501;i++)a[i]+=weiyi;}//还原A矩阵return (v2+weiyi);}double fanmifayunsuan(double *a,double weiyi){int i,k;double b=0.16;double c=-0.064;double ee,w,v1,v2,mm,sum;double u[501];double y[501];double C[5][501];void yasuo(double *A,double (*C)[501]);void LUfenjie(double (*C)[501]);void qiuU(double (*C)[501],double *y,double *u);//把A阵压缩到C阵中for(i=0;i<501;i++)u[i]=1;//给u赋初值if (weiyi!=0){for (i=0;i<501;i++)a[i]-=weiyi;}yasuo(a,C);LUfenjie(C);ee=1;k=0; //使得初始计算时进入循环语句while(ee>1e-12){mm=0;for(i=0;i<501;i++){mm=mm+u[i]*u[i];}w=sqrt(mm);for(i=0;i<501;i++){y[i]=u[i]/w;}qiuU(C,y,u);sum=0;for(i=0;i<501;i++){sum+=y[i]*u[i];}v1=v2;v2=sum;//去除特殊情况,减少漏洞if(k==0){k++;}else{ee=fabs(1/v2-1/v1)/fabs(1/v2);}}if (weiyi!=0){for (i=0;i<501;i++)a[i]+=weiyi;}//还原A矩阵return (1/v2+weiyi);}void yasuo(double *A,double (*C)[501]){double b=0.16;double c=-0.064;int i;for(i=0;i<501;i++){C[0][i]=c;C[1][i]=b;C[2][i]=A[i];C[3][i]=b;C[4][i]=c;}}void LUfenjie(double (*C)[501]){int k,t,j;int r=2,s=2;double sum;int minn(int ,int );int maxx(int ,int );for(k=0;k<501;k++){for(j=k;j<=minn(k+s,501-1);j++){if(k==0)sum=0;else{sum=0;for(t=maxx(k-r,j-s);t<k;t++){sum=sum+C[k-t+s][t]*C[t-j+s][j];}}C[k-j+s][j]=C[k-j+s][j]-sum;}for(j=k+1;j<=minn(k+r,501-1);j++){if(k<501-1){if(k==0)sum=0;else{sum=0;for(t=maxx(j-r,k-s);t<k;t++){sum=sum+C[j-t+s][t]*C[t-k+s][k];}}C[j-k+s][k]=(C[j-k+s][k]-sum)/C[s][k];}}}}void qiuU(double (*C)[501],double *y,double *u){int i,t;double b[501];double sum;int r=2,s=2;int minn(int ,int );int maxx(int ,int );for(i=0;i<501;i++){b[i]=y[i];}for(i=1;i<501;i++){sum=0;for(t=maxx(0,i-r);t<i;t++){sum=sum+C[i-t+s][t]*b[t];}b[i]=b[i]-sum;}u[500]=b[500]/C[s][500];for(i=501-2;i>=0;i--){sum=0;for(t=i+1;t<=minn(i+s,500);t++){sum=sum+C[i-t+s][t]*u[t];}u[i]=(b[i]-sum)/C[s][i];}}int minn(int x,int y){int min;if(x>y)min=y;elsemin=x;return min;}int maxx(int b,int c){int max;if(b>c){if(b>0)max=b;elsemax=0;}else{if(c>0)max=c;elsemax=0;}return max;}三、特征值以及的值λ[1]=-1.070011361502e+001 λ[501]=9.724634098777e+000λ[s]=-5.557910794230e-003最接近u[1]的特征值为-1.018293403315e+001最接近u[2]的特征值为-9.585707425068e+000最接近u[3]的特征值为-9.172672423928e+000最接近u[4]的特征值为-8.652284007898e+000最接近u[5]的特征值为-8.0934********e+000最接近u[6]的特征值为-7.659405407692e+000最接近u[7]的特征值为-7.119684648691e+000最接近u[8]的特征值为-6.611764339397e+000最接近u[9]的特征值为-6.0661********e+000最接近u[10]的特征值为-5.585101052628e+000最接近u[11]的特征值为-5.114083529812e+000最接近u[12]的特征值为-4.578872176865e+000最接近u[13]的特征值为-4.096470926260e+000最接近u[14]的特征值为-3.554211215751e+000最接近u[15]的特征值为-3.0410********e+000最接近u[16]的特征值为-2.533970311130e+000最接近u[17]的特征值为-2.003230769563e+000最接近u[18]的特征值为-1.503557611227e+000最接近u[19]的特征值为-9.935586060075e-001最接近u[20]的特征值为-4.870426738850e-001最接近u[21]的特征值为2.231736249575e-002最接近u[22]的特征值为5.324174742069e-001最接近u[23]的特征值为1.052898962693e+000最接近u[24]的特征值为1.589445881881e+000最接近u[25]的特征值为2.060330460274e+000最接近u[26]的特征值为2.558075597073e+000最接近u[27]的特征值为3.080240509307e+000最接近u[28]的特征值为3.613620867692e+000最接近u[29]的特征值为4.0913********e+000最接近u[30]的特征值为4.603035378279e+000最接近u[31]的特征值为5.132924283898e+000最接近u[32]的特征值为5.594906348083e+000最接近u[33]的特征值为6.080933857027e+000最接近u[34]的特征值为6.680354092112e+000最接近u[35]的特征值为7.293877448127e+000最接近u[36]的特征值为7.717111714236e+000最接近u[37]的特征值为8.225220014050e+000最接近u[38]的特征值为8.648666065193e+000最接近u[39]的特征值为9.254200344575e+000cond(A)2=1.925204273902e+003 Det(A)=2.772786141752e+118四、现象讨论在大作业的程序设计过程当中,初始向量的赋值我顺其自然的设为第一个分量为1,其它分量为0的向量,计算结果与参考答案存在很大差别,计算结果对比如下图2所示(左侧为正确结果,右侧为错误结果),导致了我花了很多的时间去检查程序算法。

计算实习报告三篇

计算实习报告三篇

If you think that you have forgotten a person, then you are not so stupid as to mention her with forgetting.悉心整理助您一臂(页眉可删)计算实习报告三篇计算实习报告篇1一、实习目的通过理论联系实际,巩固所学的知识,提高处理实际问题的能力,为顺利毕业进行做好充分的准备,并为自己能顺利与社会环境接轨做准备。

通过这次实习,使我们进一步理解和领会所学的基本理论,了解计算机技术和信息管理技术的发展及应用,较为系统地掌握计算机应用技能和信息管理技能,把所学知识与解决实际问题相联系,能够利用计算机处理工作中的各种信息,培养我们发现问题、分析问题和解决问题的能力,从而提高我们从事实际工作的能力。

通过理论联系实际,巩固所学的知识,提高处理实际问题的能力,了解设计专题的主要内容,使学生能够了解社会、学校的需要,在实习单位领导的帮助,对自己今后所从事的事业有一个实习了解的过程。

为毕业设计的顺利进行做好充分的准备,并为自己能顺利与社会环境接轨做准备。

二、实习意义通过实习,使学生在社会实践中接触与本专业相关的实际工作,增强感性认识,培养和锻炼学生综合运用所学的基础理论、基本技能和专业知识,去独立分析和解决实际问题的能力,把理论和实践结合起来,提高实践动手能力,为学生毕业后走上工作岗位打下一定的基础;同时可以检验教学效果,为进一步提高教育教学质量,培养合格人才积累经验。

计算机是一门对实践要求较高的学科,通过专业实习,使学生能熟悉有关计算机专业的各个领域,使学生毕业后能胜任与本专业相关的工作。

大学四年学习了很多,经历了很多,得到的是学习能力、处事能力和一些专业知识。

可面对社会,我们经验太少,思想单纯!毕业实习,给了我们一个了解社会,增加经验,熟悉工作单位的机会。

锻炼自己的动手能力,将学习的理论知识运用于实践当中,反过来还能检验书本上理论的正确性,有利于融会贯通。

数值分析实验报告

数值分析实验报告

数值分析实验报告一、实验目的数值分析是一门研究用计算机求解数学问题的数值方法及其理论的学科。

本次实验的目的在于通过实际操作和编程实现,深入理解和掌握数值分析中的常见算法,提高运用数值方法解决实际问题的能力,并对算法的精度、稳定性和效率进行分析和比较。

二、实验环境本次实验使用的编程语言为 Python,使用的开发工具为 PyCharm。

实验所依赖的主要库包括 NumPy、Matplotlib 等。

三、实验内容(一)函数逼近与插值1、拉格朗日插值法通过给定的离散数据点,构建拉格朗日插值多项式,对未知点进行函数值的估计。

2、牛顿插值法与拉格朗日插值法类似,但采用了不同的形式和计算方式。

(二)数值积分1、梯形公式将积分区间划分为若干个梯形,通过计算梯形面积之和来近似积分值。

2、辛普森公式基于抛物线拟合的方法,提高积分近似的精度。

(三)线性方程组求解1、高斯消元法通过逐行消元将线性方程组化为上三角形式,然后回代求解。

2、 LU 分解法将系数矩阵分解为下三角矩阵 L 和上三角矩阵 U,然后通过两次前代和回代求解。

(四)非线性方程求解1、二分法通过不断将区间一分为二,逐步缩小根所在的区间,直到满足精度要求。

2、牛顿迭代法利用函数的切线来逼近根,通过迭代逐步收敛到根的近似值。

四、实验步骤(一)函数逼近与插值1、拉格朗日插值法定义计算拉格朗日基函数的函数。

根据给定的数据点和待求点,计算插值多项式的值。

输出插值结果,并与真实值进行比较。

2、牛顿插值法计算差商表。

构建牛顿插值多项式。

进行插值计算和结果分析。

(二)数值积分1、梯形公式定义积分区间和被积函数。

按照梯形公式计算积分近似值。

分析误差。

2、辛普森公式同样定义积分区间和被积函数。

运用辛普森公式计算积分近似值。

比较与梯形公式的精度差异。

(三)线性方程组求解1、高斯消元法输入系数矩阵和右端项向量。

进行消元操作。

回代求解方程。

输出解向量。

2、 LU 分解法对系数矩阵进行 LU 分解。

《数值分析》课程实验报告范文

《数值分析》课程实验报告范文

《数值分析》课程实验报告范文《数值分析》课程实验报告姓名:学号:学院:机电学院日期:2022年某月某日目录实验一函数插值方法1实验二函数逼近与曲线拟合5实验三数值积分与数值微分7实验四线方程组的直接解法9实验五解线性方程组的迭代法15实验六非线性方程求根19实验七矩阵特征值问题计算21实验八常微分方程初值问题数值解法24实验一函数插值方法一、问题提出对于给定的一元函数的n+1个节点值。

试用Lagrange公式求其插值多项式或分段二次Lagrange插值多项式。

实验二函数逼近与曲线拟合一、问题提出从随机的数据中找出其规律性,给出其近似表达式的问题,在生产实践和科学实验中大量存在,通常利用数据的最小二乘法求得拟合曲线。

在某冶炼过程中,根据统计数据的含碳量与时间关系,试求含碳量与时间t的拟合曲线。

t(分)051015202530354045505501.272.162.863.443.874.154.374.51 4.584.024.64二、要求1、用最小二乘法进行曲线拟合;2、近似解析表达式为;3、打印出拟合函数,并打印出与的误差,;4、另外选取一个近似表达式,尝试拟合效果的比较;5、某绘制出曲线拟合图。

三、目的和意义1、掌握曲线拟合的最小二乘法;2、最小二乘法亦可用于解超定线代数方程组;3、探索拟合函数的选择与拟合精度间的关系四、实验步骤:第一步先写出线性最小二乘法的M文件functionc=lpoly(某,y,m)n=length(某);b=zero(1:m+1);f=zero(n,m+1); fork=1:m+1f(:,k)=某.^(k-1);enda=f'某f;b=f'某y';c=a\b;c=flipud(c);第二步在命令窗口输入:>>lpoly([0,5,10,15,20,25,30,35,40,45,50,55],[0,1.27,2.16,2.86,3.44,3.87,4.15,4.37,4.51,4.58,4.02,4.64],2)回车得到:an=-0.00240.20370.2305即所求的拟合曲线为y=-0.0024某2+0.2037某+0.2305在编辑窗口输入如下命令:>>某=[0,5,10,15,20,25,30,35,40,45,50,55];>>y=-0.0024某某.^2+0.2037某某+0.2305;>>plot(某,y)命令执行得到如下图五、实验结论分析复杂实验数据时,常采用分段曲线拟合方法。

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

航空航天大学《数值分析》计算实习报告第一大题学院:自动化科学与电气工程学院专业:控制科学与工程学生姓名:学号:教师:电话:完成日期: 2015年11月6日航空航天大学Beijing University of Aeronautics and Astronautics实习题目:第一题 设有501501⨯的实对称矩阵A ,⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡=5011A a b c b c c b c b a其中,064.0,16.0),501,,2,1(64.0)2.0sin()024.064.1(1.0-==⋅⋅⋅=--=c b i e i i a ii 。

矩阵A 的特征值为)501,,2,1(⋅⋅⋅=i i λ,并且有||min ||,501150121i i s λλλλλ≤≤=≤⋅⋅⋅≤≤1.求1λ,501λ和s λ的值。

2.求A 的与数4015011λλλμ-+=kk 最接近的特征值)39,,2,1(⋅⋅⋅=k k i λ。

3.求A 的(谱数)条件数2)A (cond 和行列式detA 。

说明:1.在所用的算法中,凡是要给出精度水平ε的,都取12-10=ε。

2.选择算法时,应使矩阵A 的所有零元素都不储存。

3.打印以下容: (1)全部源程序;(2)特征值),,39,...,2,1(,s 5011=k k i λλλλ以及A det ,)A (cond 2的值。

4.采用e 型输出实型数,并且至少显示12位有效数字。

一、算法设计方案1、求1λ,501λ和s λ的值。

由于||min ||,501150121i i s λλλλλ≤≤=≤⋅⋅⋅≤≤,可知绝对值最大特征值必为1λ和501λ其中之一,故可用幂法求出绝对值最大的特征值λ,如果λ=0,则1λ=λ,否则501λ=λ。

将矩阵A 进行一下平移:I -A A'λ= (1)对'A 用幂法求出其绝对值最大的特征值'λ,则A 的另一端点特征值1λ或501λ为'λ+λ。

s λ为按模最小特征值,||min ||5011i i s λλ≤≤=,可对A 使用反幂法求得。

2、求A 的与数4015011λλλμ-+=kk 最接近的特征值)39,...,2,1(=k k i λ。

计算1)1,2,...,50=(i i λ-k μ,其模值最小的值对应的特征值k λ与k μ最接近。

因此对A 进行平移变换:)39,,2,1k -A A k k ==(I μ (2)对k A 用反幂法求得其模最小的特征值'k λ,则k λ='k λ+k μ。

3、求A 的(谱数)条件数2)(A cond 和行列式detA 。

由矩阵A 为非奇异对称矩阵可得:||)(min max2λλ=A cond (3)其中max λ为按模最大特征值,min λ为按模最小特征值,通过第一问我们求得的λ和s λ可以很容易求得A 的条件数。

在进行反幂法求解时,要对A 进行LU 分解得到。

因L 为单位下三角阵,行列式为1,U 为上三角阵,行列式为主对角线乘积,所以A 的行列式等于U 的行列式,为U 的主对角线的乘积。

二、 算法实现1、矩阵存储原矩阵A 为一个上、下半带宽都为2的501×501的带状矩阵,由于矩阵中的0元素太多,如果分配一个501×501的空间保存矩阵的话会浪费很多空间。

因此,为了节省存储量,A 的带外元素不给存储,值存储带元素,如下C 矩阵所示:⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡=000000501500499321cccc b b b b b a a a a a a b bb b bc c c c C(4) C 是一个5×501的矩阵,相比A 大大节省了存储空间,在数组C 中检索矩阵A 的带元素ij a 的方法是:j 1,s j -i c a A ++=中的元素的带内元素C ij (5)2、幂法幂法迭代公式如下:⎪⎪⎪⎩⎪⎪⎪⎨⎧====∈-------k T k kk k k k k u y Ay u u y u u 111k 11211k n0/R βηη任取非零向量 (6)其中λβ=∞→k k lim ,不断迭代当εβββ≤--||/||1k k k 时即可认为其满足精度要求,令k βλ=。

在程序中计算1u -=k k Ay 时,根据A 矩阵的特点,简化如下:⎪⎪⎪⎪⎩⎪⎪⎪⎪⎨⎧⋅⋅⋅=+++++-+-=++=+++=+++=++=)499,,3()2()1()()()1()2()()501()500()499()501()501()500()499()498()500()4()3()2()2(C )1()2()3()2()1(C )1(]][3[]501][3[]500][3[]2][3[]1][3[i i cy i by i y i C i by i cy i u y C by cy u by y C by cy u cy by y by u cy by y u i (7)3、反幂法反幂法迭代公式如下:⎪⎪⎪⎩⎪⎪⎪⎨⎧====∈-------kT k k k k k k k u y y Au u y u u 111k 11211k n0/R βηη任取非零向量 (8)当k 足够大时,kβλ1s ≈。

在求解1-=k k y Au 时,可先对A 进行Doolittle 分解,由于A 是带状结构,所以分解出的L 、U 也是带状结构,利用C 矩阵进行Doolittle 分解并求向量k u 的算法如下: (1)作分解A=LU对于n ,,2,1k =执行:)],min(,,2,1[/)(:)],min(,,1,[:,11),,1max(,1,1,1,11),,1max(,1,1,1,1n r k k k i c ccc c n s k k k j ccc c ks k s k r i t k s k t t s t i k s k i k s k i k s j r k t js j t t s t k j s j k j s j k +++=-=++=-=+---=++-++-++-++----=++-++-++-++-∑∑ (9)由于C 语言中数组下标是从0开始的,所以在程序中矩阵元素c 的下标都减1。

(2)求解y Ux b Ly ==,(数组b 先是存放原方程组右端向量,后来存放中间向量y,在程序中b 和y 都保存在数组y[501]中。

))1,,2,1(/)(/),,3,2(,1),min(1,1,11),1max(,1 --=-===-=+++=++-+--=++-∑∑n n i c x cb xc b x n i b cb b i s tn s i i t t s t i i i ns n n i r i t tt s t i i i (10)求出k u 后,其他部分与幂法求解相同。

三、结果分析实验表明,本程序中,初始向量[]Tu u u u 501211501,,, =⨯对结果影响较大,合适的初始向量对得到正确的收敛结果比较重要,如表1是不同初始向量的情况下的得到的部分结果。

(实验结果截图见附录)表1 不同初始向量对应的1λ和501λ及其迭代次数由表1可以得到如下结论:1. 不同的初始向量对本程序的1λ影响大,对501λ没有影响,都能保证收敛到正确值。

2. 初始向量中必须保证501462~u u 中至少有一个为1才能保证1λ收敛到正确值。

3. 初始向量非零值的多少和大小对迭代次数并没有明显影响。

4.为解决初始向量对程序的影响,可以先对A做平移变换再求1四、实验程序#include <stdio.h>#include <math.h>#include <stdlib.h>static double b=0.16,c=-0.064;#define Precision 1e-12void copy(double b[501],double y[501]);double dianji(double x[],double y[]);//计算两个向量积void InitMatrix(double *p);//Init Matrix Adouble NeiJi(double a[],double b[]);//get 2数void get_y(double *y,double *u);//get yvoid get_u(double *u,double *y,double *a);//get uvoid Initu(double *p);//初始化初始向量udouble Get_Fabs_Eigenvalue(double *a,double *u,int *iterations);//循环迭代得到绝对值最大特征值void A_sub_minI(double *a,double min);//A-min*Ivoid InitC(double C[5][501]);//初始化数组Cvoid DoolittleC(double C[5][501],int n,int s,int r);//进行Doolittle分解int min(int a,int b);//取返回a,b最小值int max(int a,int b);//求最大值void Doolittle_getx(double C[5][501],double y[501],double u[501],int n,int s,int r);double Get_min_Eigenvalue(double C[5][501],double *u,int n,int s,int r,int *iterations);double detA(double C[5][501]);//求A的行列式struct FinalValue{double min;//特征值最小值double max;//特征值最大值double abs_min;//模最小特征值double abs_max;//模最大特征值double detA;//A的行列式double cond2;//A的条件数int min_iterations;//最小值迭代次数int max_iterations;//最大值迭代次数int abs_min_iterations;//求绝对值最小的迭代次数};int main(){FinalValue main_num= {0,0,0,0};double temp;//两值交换中间变量int temp1;double u[501]={0};//,y[501];//为u0赋初值double Norm_u=0;//数double C[5][501] = {0};InitC(C);//初始化Cint i=0,*iterations;Initu(u);iterations = &main_num.min_iterations;main_num.min = Get_Fabs_Eigenvalue(C[2],u,iterations);//将求得的绝对值最大特征值放到min变量中main_num.abs_max = main_num.min;A_sub_minI(C[2],main_num.min);for(i=0;i<501;i++)u[i]=i+1;iterations = &main_num.max_iterations;main_num.max = Get_Fabs_Eigenvalue(C[2],u,iterations);//将求得的绝对值最大特征值放到min变量中main_num.max += main_num.min;if(main_num.min>main_num.max){temp = main_num.min;main_num.min = main_num.max;main_num.max = temp;temp1 = main_num.min_iterations;main_num.min_iterations = main_num.max_iterations;main_num.max_iterations = temp1;}printf("最小特征值为:%.12e,迭代次数为:%d\n",main_num.min,main_num.min_iterations);printf("最大特征值为:%.12e,迭代次数为:%d\n",main_num.max,main_num.max_iterations);/**********************************//*以下利用反幂法求解模最小的特征值*//**********************************/InitC(C);//初始化CInitu(u);DoolittleC(C,501,2,2);main_num.detA = detA(C);iterations = &main_num.abs_min_iterations;main_num.abs_min = Get_min_Eigenvalue(C,u,501,2,2,iterations);main_num.cond2 = fabs(main_num.abs_max/main_num.abs_min);printf("绝对值最小特征值为:%.12e,迭代次数为:%d\n",main_num.abs_min,main_num.abs_min_iterations);printf("A的行列式为:%.12e;",main_num.detA);printf("A的条件数cond(A)2为:%.12e\n",main_num.cond2);/**********************************************//*以下利用反幂法求与数列Uk中元素最相近的特征值*//**********************************************/double Uk[39];//保存Uk的值并保存与Uk[i]最接近的λdouble B[39];int diedai;iterations = &diedai;for( i=1;i<=39;i++){Uk[i-1] = main_num.min + i*(main_num.max - main_num.min)/40;InitC(C);Initu(u);for(int j =0;j<501;j++)C[2][j] -= Uk[i-1];DoolittleC(C,501,2,2);B[i-1] = Get_min_Eigenvalue(C,u,501,2,2,iterations);B[i-1] +=Uk[i-1];printf("μ%-2d=%-3.12e,与μ%-2d最相近的特征值λ=%-3.12e\n",i,Uk[i-1],i,B[i-1]);}return 0;}double detA(double C[5][501]){int i =0;double e =1;for(i =0;i<501;i++)e*=C[2][i];return e;}void InitMatrix(double *p){for(int i=1;i<=501;i++){*p=(1.64-0.024*i)*sin(0.2*i)-0.64*exp(0.1/i);p++;}}void Initu(double *p){for(int i=1;i<=501;i++){if(i<=500)*p=0;else *p=1;p++;}// *p=1;}void A_sub_minI(double *a,double min) {for(int i=1;i<=501;i++){*a=*a-min;a++;}}double NeiJi(double *a,double *b) {double e=0.0;for(int i=0;i<501;i++){e=e+(*a)*(*b);a++;b++;}return e;}void get_y(double *y,double *u){ double Norm_u=sqrt(NeiJi(u,u));for(int i=0;i<501;i++){*y = *u/Norm_u;y++;u++;}}void get_u(double *u,double *y,double *a){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];for(int i=2;i<499;i++)u[i]=c*y[i-2]+b*y[i-1]+a[i]*y[i]+b*y[i+1]+c*y[i+2];u[499]=c*y[497]+b*y[498]+a[499]*y[499]+b*y[500];u[500]=c*y[498]+b*y[499]+a[500]*y[500];}void InitC(double C[5][501]){int i;for(i = 2;i < 501;i++){C[0][i] = -0.064;C[4][i-2] = -0.064;}for(i = 1;i < 501;i++){C[1][i] = 0.16;C[3][i-1] = 0.16;}for(i = 1;i <= 501;i++){C[2][i-1] = (1.64-0.024*i)*sin(0.2*i)-0.64*exp(0.1/i);}}double Get_Fabs_Eigenvalue(double *a,double *u,int *iterations) {double y[501],B_k0=0,B_k1=0;double wucha;int i=0;while(1){++i;get_y(y,u);get_u(u,y,a);B_k1=NeiJi(y,u);//是否判断B_K1是否为0?wucha=(fabs(B_k1-B_k0))/(fabs(B_k1));//get wuchaif(wucha<=Precision)break;else if(i>10000){printf("迭代次数超长,请更改初始向量\n");break;}else B_k0 = B_k1;}*iterations = i;return B_k1;}double Get_min_Eigenvalue(double C[5][501],double *u,int n,int s,int r,int *iterations){double y[501] = {0},B_k0=0,B_k1=0;double b[501] = {0};//储存y值的中间向量double wucha;int i=0;while(1){++i;get_y(y,u);copy(b,y);//保护y向量Doolittle_getx(C,y,u,501,2,2);copy(y,b);B_k1=NeiJi(y,u);//是否判断B_K1是否为0?wucha=(fabs(1/B_k1-1/B_k0))/(1/fabs(B_k1));//get wucha// wucha=(fabs(B_k1-B_k0))/(fabs(B_k1));//get wuchaif(wucha<=Precision)break;else if(i>10000){printf("迭代次数超长,请更改初始向量\n");break;}else B_k0 = B_k1;}*iterations = i;return (1/B_k1);}void Doolittle_getx(double C[5][501],double y[501],double u[501],int n,int s,int r){int i,t;double e;for(i = 2;i <= n;i++){e = 0;for(t = max(1,i-r);t <= i-1;t++){e+=C[i-t+s+1-1][t-1]*y[t-1];}y[i-1] -= e;}u[n-1] = y[n-1]/C[s+1-1][n-1];for(i = n-1;i>=1;i--){e = 0;for(t = i+1;t<=min(i+s,n);t++)e += C[i-t+s+1-1][t-1]*u[t-1];u[i-1] = (y[i-1] - e)/C[s+1-1][i-1];}}void copy(double b[501],double y[501]){for(int i=0;i<501;i++){b[i] = y[i];}}void DoolittleC(double C[5][501],int n,int s,int r){int k,j,i,t;double e;for(k = 1;k <= n;k++){for(j = k;j <= min(k+s,n);j++){ e = 0;for(t = max(max(1,k-r),j-s);t <= k-1;t++){e = e + C[k-t+s+1-1][t-1]*C[t-j+s+1-1][j-1];}C[k-j+s+1-1][j-1] = C[k-j+s+1-1][j-1] - e;}if(k == n)break;for(i = k+1;i <= min(k+r,n);i++){ e = 0;for(t = max(max(1,i-r),k-s);t <= k-1;t++){e = e + C[i-t+s+1-1][t-1]*C[t-k+s+1-1][k-1];}C[i-k+s+1-1][k-1] = (C[i-k+s+1-1][k-1] - e)/C[s+1-1][k-1];}}}int min(int a,int b){if(a<b)return a;else return b;}int max(int a,int b){if(a>=b)return a;elsereturn b;}附录:部分实验程序截图1、[]Tu 1,,1,15011 =⨯2、[]Tu 501,,2,11501 =⨯3、[]Tu 0,,0,11501 =⨯4、[]Tu 1,0,,01501 =⨯5、[]Tu 0,0,1,11501 ,=⨯6、146121====u u u7、140021====u u u8、另14621===u u ,其余都为09、1462 u ,其余为010、1501462===u u ,其余为011、1481 u ,其余为0。

相关文档
最新文档