数值方法实验报告.doc

2009-2010学年 第一学期数值方法实验报告 学 院： 石油工程学院专业班级： 指导教师： 李梦霞学生姓名： 实验一1.题目 高斯列主元素消去法：用Gauss列主元消去法解线性方程组，其中，A= 2.理论分析及算法描述（1）理论分析： 由高斯消去法，当=0时，消去过程则中断。此外，即使0，但当绝对值较小时，也会因用它作除数而出现舍入误差急剧增长的系数，从而导致计算结果不准确。 列主元高斯消去法只需在高斯消去法消去过程的第k步消去计算前，插入搜索绝对值最大元和交换过程（交换增广矩阵的行）的处理。（2）Guass 列主元消去法算法描述：Step1:k=1；Step2:按列选主元|=|，；Step3:if|=0，停止计算；Step4:=K，则不换行。Else，交换行与K行，（j=k,.N),;Step5:计算乘数(i=k+1,.n);Step6:消元计算=-lik*=-lik*(i,j=k+1,.n);Step7:k=k+1，If k=n-1转Step2；Step8:回代求解（k=n,n-1,.1）。3.源程序#include #include #define N 6main() int i,j,k; float p,m; float aNN,bN,lN; FILE *fp1=fopen(d:date1.txt,rw); FILE *fp2=fopen(d:date2.txt,rw); for(i=0;iN;i+) fscanf(fp2,%f,&bi); printf(%8.3fn,bi); fclose(fp2); for(i=0;iN;i+) for(j=0;jN;j+) fscanf(fp1,%f,&aij); printf(%8.3f ,aij); printf(n); fclose(fp1); for(j=0;jN;j+) for(i=j+1;iabs(ajj) for(k=0;kN;k+) p=aik; m=bi; aik=ajk; bi=bj; ajk=p; bj=m; if(ajj=0) printf(can not do it); break; for(j=0;jN;j+) for(i=j+1;iN;i+) li=aij/ajj; for(k=j+1;k=0;i-) m=0.000; for(k=i;kN;k+) m=m+aik*bk; bi=(bi-m)/aii; for(i=0;iN;i+) printf(X%d解为：,i+1);printf(%8.3fn,bi); 4.运行结果 实验二1.题目 用Cholesky分解，解下列方程组,其中： 2.理论分析及算法描述(1)理论分析：在工程计算中，常遇到所求解的线性方程组的系数矩阵具有对称正定性，矩阵的这一特性使它的三角分解具有更简单的形式。设矩阵A为对称正定矩阵，即A的各阶顺序主子式均大于零，则由定理可知存在一个单位下三角矩阵和一个上三角矩阵U使得 A=U 成立。定理2.3 如果A为对称正定矩阵，则存在唯一的对角线元素全为正数的下三角形矩阵L，使得 A=L 成立。当矩阵A完成Cholesky分解后，求解方程组Ax=b就转化为求解方程组 Ly=b , x=y 。(2)算法描述：Step1：Step2：Step3：3.源程序Chosy分解#include #include #define N 8main() int i,j,k; float m; float aNN,bN; FILE *fp1=fopen(d:date3.txt,rw); FILE *fp2=fopen(d:date4.txt,rw); for(i=0;iN;i+) for(j=0;jN;j+) fscanf(fp1,%f,&aij);printf(%4.1f ,aij); printf(n); fclose(fp1); for(i=0;iN;i+) fscanf(fp2,%f,&bi); printf(%4.1fn,bi); fclose(fp2); a00=(float)sqrt(a00); for(i=1;iN;i+) for(j=0;ji;j+) m=0.0; for(k=0;kj;k+)m=aik*ajk; aij=(aij-m)/ajj; m=0.0; for(k=0;ki;k+)m=aik*aik;aii=(float)sqrt(aii-m); b0=b0/a00; for(i=1;iN;i+) m=0.0; for(k=0;k=0;i-) m=0.0; for(k=i+1;kN;k+) m=aik*bk; bi=(bi-m)/aii; for(i=0;i|q| 则q=P (3)=+PStep5:if|q|e转Step2 否则转Step6Step6：输出，k。3.源程序#include #include #define N 10main() int i,j,k; float p,m,w,t,q=0.0; double e; float aNN,bN,xN=1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0; FILE *fp1=fopen(d:date5.txt,rw); FILE *fp2=fopen(d:date6.txt,rw); printf(请输入精度要求); scanf(%f,&e); for(i=0;iN;i+) for(j=0;jN;j+) fscanf(fp1,%f,&aij);printf(%4.1f ,aij); printf(n); fclose(fp1); for(i=0;iN;i+) fscanf(fp2,%f,&bi); printf(%4.1fn,bi); fclose(fp2); printf(请输入松弛因子); scanf(%f,&w); q=0.0; do k+; m=0.0; t=0.0; for(i=0;iN;i+) for(j=0;ji-1;j+) m=m+aij*xj; for(j=i-1;jfabs(q) q=p; xi=xi+q; while(fabs(q)e&k12); for(i=0;iN;i+) printf(%4.1fn,xi); printf(迭代次数 %d,k);4.运行结果 实验四1.题目 已知当时，对应的函数值=1，4，7，8，6，求四次牛顿插值多项式。2.理论分析及算法描述(1)理论分析： 其中 可以看出是x的n次多项式，并且当时，有，于是 (i=0,1,.,n)亦即满足插足条件，因此它就是所要构造的Newton插值公式，称为Newton插值公式的余项。(2)算法描述：Step1：输入Step2： （k=i,.,N-1）Step3：if i=N 转step4Step4：输出3.源程序#include #include #define N 5main() float aN,bN; int i,k; FILE *fp1=fopen(d:date7.txt,rw); FILE *fp2=fopen(d:date8.txt,rw); printf(x的值为n); for(i=0;iN;i+) fscanf(fp1,%f,&ai);printf(%4.1f ,ai); printf(nn); fclose(fp1); printf(y的值为n); for(i=0;iN;i+) fscanf(fp2,%f,&bi); printf(%4.1f,bi); printf(nn); fclose(fp2); for(i=1;iN;i+) for(k=i;k,i=0,转Step7Step4:若=0，i=0,转Step7Step5:计算，Step6：若|-|e或|e,则输出 否则=，=，转Step2Step7：当i=0时，输出失败信息。3.源程序#include #include #define N 10main() int i,k; float a,b,c,d,e,m,t,x,p; printf(请输入相应的已知系数); scanf(%f %f %f %f %f,&a,&b,&c,&d,&e); printf(请输入精度要求); scanf(%f,&p); x=5.0; for(k=1;kN;k+) m=0.0; t=0.0; m=a*x*x*x*x+b*x*x*x+c*x*x