西交计算方法A上机大作业
计算方法上机报告_董晓壮

计算方法A上机报告学院(系):电气工程学院学生姓名:陶然学号:授课老师:完成日期:2019年12月03日西安交通大学Xi'an Jiaotong University目录1 QR分解法求解线性方程组 (2)1.1 算法原理 (2)1.1.1 基于吉文斯变换的QR分解 (2)1.1.2 基于豪斯霍尔德变换的QR分解 (3)1.2 程序流程图 (4)1.2.1 基于吉文斯变换的QR分解流程图 (4)1.2.2 基于豪斯霍尔德变换的QR分解流程图 (5)1.3 程序使用说明 (5)1.3.1 基于吉文斯变换的QR分解程序说明 (5)1.3.2 基于豪斯霍尔德变换的QR分解程序说明 (7)1.4 算例计算结果 (8)2 共轭梯度法求解线性方程组 (10)2.1 算法原理 (10)2.2 程序流程图 (10)2.3 程序使用说明 (11)2.4 算例计算结果 (12)3 三次样条插值 (14)3.1 算法原理 (14)3.2 程序流程图 (16)3.3 程序使用说明 (17)3.4 算例计算结果 (19)4 龙贝格积分 (21)4.1 算法原理 (21)4.2 程序流程图 (22)4.3 程序使用说明 (23)4.4 算例计算结果 (24)结论 (26)1 QR 分解法求解线性方程组1.1 算法原理矩阵的QR 分解是指,可以将矩阵A 分解为一个正交矩阵Q 和一个上三角矩阵R 的乘积,实际中,QR 分解经常被用来解决线性最小二乘问题,分解情况如图1.1所示。
=⨯图1.1 QR 分解示意图本次上机学习主要进行了两个最基本的正交变换—吉文斯(Givens )变换和豪斯霍尔德(Householder )变换,并由此导出矩阵的QR 分解以及求解线性方程组的的方法。
1.1.1 基于吉文斯变换的QR 分解吉文斯矩阵也称初等旋转阵,如式(1.1)所示,它把n 阶单位矩阵I 的第,i j 行的对角元改为c ,将第i 行第j 列的元素改为s ,第j 行第i 列的元素改为s −后形成的矩阵。
西安交通大学计算方法B大作业

计算方法上机报告姓名:学号:班级:目录题目一------------------------------------------------------------------------------------------ - 4 -1.1题目内容 ---------------------------------------------------------------------------- - 4 -1.2算法思想 ---------------------------------------------------------------------------- - 4 -1.3Matlab源程序----------------------------------------------------------------------- - 5 -1.4计算结果及总结 ------------------------------------------------------------------- - 5 - 题目二------------------------------------------------------------------------------------------ - 7 -2.1题目内容 ---------------------------------------------------------------------------- - 7 -2.2算法思想 ---------------------------------------------------------------------------- - 7 -2.3 Matlab源程序---------------------------------------------------------------------- - 8 -2.4计算结果及总结 ------------------------------------------------------------------- - 9 - 题目三----------------------------------------------------------------------------------------- - 11 -3.1题目内容 --------------------------------------------------------------------------- - 11 -3.2算法思想 --------------------------------------------------------------------------- - 11 -3.3Matlab源程序---------------------------------------------------------------------- - 13 -3.4计算结果及总结 ------------------------------------------------------------------ - 14 - 题目四----------------------------------------------------------------------------------------- - 15 -4.1题目内容 --------------------------------------------------------------------------- - 15 -4.2算法思想 --------------------------------------------------------------------------- - 15 -4.3Matlab源程序---------------------------------------------------------------------- - 15 -4.4计算结果及总结 ------------------------------------------------------------------ - 16 - 题目五----------------------------------------------------------------------------------------- - 18 -5.1题目内容 --------------------------------------------------------------------------- - 18 -5.2算法思想 --------------------------------------------------------------------------- - 18 -5.3 Matlab源程序--------------------------------------------------------------------- - 18 -5.3.1非压缩带状对角方程组------------------------------------------------- - 18 -5.3.2压缩带状对角方程组---------------------------------------------------- - 20 -5.4实验结果及分析 ------------------------------------------------------------------ - 22 -5.4.1Matlab运行结果 ---------------------------------------------------------- - 22 -5.4.2总结分析------------------------------------------------------------------- - 24 -5.5本专业算例 ------------------------------------------------------------------------ - 24 - 学习感悟-------------------------------------------------------------------------------------- - 27 -题目一1.1题目内容计算以下和式:0142111681848586n n S n n n n ∞=⎛⎫=--- ⎪++++⎝⎭∑,要求: (1)若保留11个有效数字,给出计算结果,并评价计算的算法; (2)若要保留30个有效数字,则又将如何进行计算。
西安电子科技大学出版社计算方法上机答案

西安电子科技大学出版社《计算方法》任传祥等编著第九章计算方法上机参考答案实验一,算法一#include <stdio.h>#include <math.h>double I0=log(6)/log(5),I1;int n=1;main (){while(1){I1=1.0/(n)-I0*5.0;printf("%d %lf\n", n,I1);if(n>=20)break;elseI0=I1;n++;}}实验一,算法二#include <stdio.h>#include <math.h>double I0=(1/105.0+1/126.0)/2,I1;int n=20;main (){printf("%d %lf\n", n,I0);while(1){I1=1.0/(5.0*n)-I0/5.0;printf("%d %lf\n", n-1,I1);if(n<2)break;elseI0=I1;n--;}}实验二,二分法#include <stdio.h>#include <math.h>#define esp 1e-3double f(double x);main (){double a=1,b=2,x;while(fabs(b-a)>esp){x=(a+b)/2;printf("x=%lf\n",x);if(f(x)==0)break;elseif(f(x)*f(a)<0)b=x;elsea=x;}}double f(double x){return pow(x,3)-x-1;}实验二,牛顿迭代法#include<stdio.h>#include<math.h>double f(double x);double f1(double x);#define esp 1e-3void main(){double x0 = 1.5, x1;x1 = x0 - f(x0) / f1(x0);printf("x=%lf\n", x1);x0 = x1;x1 = x0 - f(x0) / f1(x0);printf("x=%lf\n", x1);while (fabs(x1 - x0)>esp){x0 = x1;x1 = x0 - f(x0) / f1(x0);printf("x=%lf\n", x1);} }double f(double x){return pow(x, 3) - x - 1;} double f1(double x){return 3 * x*x - 1;}弦割法#include<stdio.h>#include<math.h>double f(double x);#define esp 1e-3void main(){double x0 = 1.5, x1=2.0,x2;do{ x2=x1 - (x1-x0)*f(x1) /(f(x1)-f(x0));x0=x1;x1=x2;printf("x=%lf\n", x1);}while (fabs(x1 - x0)>esp);{printf("x=%lf\n", x1);}}double f(double x){return pow(x, 3) - x - 1;}实验3#include <stdio.h>/*列主元高斯消去法*/#include <math.h>float x[3],temp,max;float A[3][4]={10,-2,-1,3,-2,10,-1,15,-1,-2,5,10},c[3][4]={10,-2,-1,3,-2,10,-1,15,-1,-2,5,10}; int n=3,i,k,j,m;void main(){for(i=0;i<n;i++){max=A[i][i];k=i;for(j=j+1;j<n;j++){{max=fabs(A[j][i]);k=j;}}if(k!=i){for(j=i+1;j<=n;j++){temp=A[i][j];A[i][j]=A[k][j];A[k][j]=temp;}}for(j=i+1;j<n;j++)for(m=i+1;m<=n;m++){c[j][m]=c[j][m]+(-c[j][i]/c[i][i])*c[i][m];}}for(i=n-1;i>=0;i--){temp=0.0;for(j=n-1;j>=i+1;j--)temp=temp+c[i][j]*x[j];x[i]=(c[i][n]-temp)/c[i][i];}printf("x[1]=%f\nx[2]=%f\nx[3]=%f\n",x[0],x[1],x[2]);实验四,拉格朗日插值#include<stdio.h>int n=5,i,j;double l,L=0,X=0.5;main(){double x[5]={0.4,0.55,0.65,0.8,0.9};doubley[5]={0.41075,0.57815,0.69675,0.88811,1.02652}; for(i=0;i<n;i++){l=y[i];for(j=0;j<n;j++){if(j!=i)l=l*(X-x[j])/(x[i]-x[j]); } L=L+l;}printf("%lf\n",L);return 0;} X=0.5 X=0.7 X=0.85牛顿插值法#include<stdio.h>#include<math.h>main(){double x[5]={0.4,0.55,0.65,0.8,0.9};doubley[5]={0.41075,0.57815,0.69675,0.88811,1.02652};int n=5,i,j;double z;printf("input z\n");scanf("%lf",&z);double a[5][5];for(i=0;i<5;i++)a[i][0]=y[i];for(i=1;i<5;i++)for(j=i;j<5;j++)a[j][i]=(a[j][i-1]-a[j-1][i-1])/(x[j]-x[j-i]);double N=a[0][0],temp=1.0;for(i=1;i<n;i++){temp=temp*(z-x[i-1]);N=N+a[i][i]*temp;}printf("N=%lf\n",N);return 0;}实验五曲线拟合#include <stdio.h>#include <math.h>float x[5]={1,2,3,4,5};float y[5]={7,11,17,27,40};float A[2][3],c[2][3];float z[2],temp,max;int i,j,k,m;int n=2;void main(){for(i=0;i<5;i++){c[0][0]=A[0][0]+=1;c[0][1]=A[0][1]+=x[i];c[0][2]=A[0][2]+=y[i];c[1][0]=A[1][0]+=x[i];c[1][1]=A[1][1]+=x[i]*x[i];c[1][2]=A[1][2]+=x[i]*y[i];}/* for(i=0;i<2;i++){printf(" %lf %lf %lf\n",A[i][0],A[i][1],A[i ][2]);}*/for(i=0;i<n;i++){max=A[i][i];k=i;for(j=j+1;j<n;j++){if(fabs(A[j][i])>max){max=fabs(A[j][i]);k=j;}} if(k!=i){for(j=i+1;j<=n;j++){temp=A[i][j];A[i][j]=A[k][j];A[k][j]=temp;}}for(j=i+1;j<n;j++)for(m=i+1;m<=n;m++){c[j][m]=c[j][m]+(-c[j][i]/c[i][i])*c[i][m];}}for(i=n-1;i>=0;i--){temp=0.0;for(j=n-1;j>=i+1;j--)temp=temp+c[i][j]*z[j];z[i]=(c[i][n]-temp)/c[i][i];}printf("a=%f\nxb=%f\n",z[0],z[1]); }实验六数值积分/*梯形*/#include<stdio.h>#include<math.h> double f(double x); main(){double x[10],y[10];double h,b=1,a=0,I;int n,i;printf("n\n");scanf("%d",&n);h=(b-a)/n;for(i=0;i<=n;i++){x[i]=a+(i*h);y[i]=f(x[i]);}I=f(a)+f(b);for(i=1;i<=n-1;i++){I=I+2*y[i];}I=(h/2)*I;printf("%lf",I);}double f(double x){double f;f=1.0/(1.0+(x*x));return(f);}/*辛普森*/#include<stdio.h>#include<math.h>double f(double x);main(){double x[30],y[30];double h,b=1,a=0,I;int n,i;printf("n\n");scanf("%d",&n);//点乘2扩展h=(b-a)/n;x[10]=1;y[10]=f(x[10]);for(i=0;i<n;i++){x[2*i]=a+(i*h);y[2*i]=f(x[2*i]);x[2*i+1]=a+(i+(1.0/2.0))*h;y[(2*i)+1]=f(x[(2*i)+1]);}I=f(a)+f(b);for(i=0;i<n;i++){I=I+4*y[(2*i)+1];}for(i=1;i<n;i++){I=I+2*y[2*i];}I=(h/6)*I;printf("%lf\n",I);}double f(double x){double f;f=1.0/(1.0+(x*x));return(f);}/*梯形*//*辛普森*/。
西安交通大学算法上机实验报告

《计算机算法设计与分析》上机实验报告姓名:班级:学号:日期:2016年12月23日算法实现题3-14 最少费用购物问题★问题描述:商店中每种商品都有标价。
例如,一朵花的价格是2元,一个花瓶的价格是5元。
为了吸引顾客,商店提供了一组优惠商品价。
优惠商品是把一种或多种商品分成一组,并降价销售。
例如,3朵花的价格不是6元而是5元。
2个花瓶加1朵花的优惠价格是10元。
试设计一个算法,计算出某一顾客所购商品应付的最少费用。
★算法设计:对于给定欲购商品的价格和数量,以及优惠价格,计算所购商品应付的最少费用。
★数据输入:由文件input.txt提供欲购商品数据。
文件的第1行中有1个整数B(0≤B≤5),表示所购商品种类数。
在接下来的B行中,每行有3个数C,K和P。
C表示商品的编码(每种商品有唯一编码),1≤C≤999;K表示购买该种商品总数,1≤K≤5;P是该种商品的正常单价(每件商品的价格),1≤P≤999。
请注意,一次最多可购买5*5=25件商品。
由文件offer.txt提供优惠商品价数据。
文件的第1行中有1个整数S(0≤S≤99),表示共有S种优惠商品组合。
接下来的S行,每行的第1个数描述优惠商品组合中商品的种类数j。
接着是j个数字对(C,K),其中C是商品编码,1≤C≤999;K表示该种商品在此组合中的数量,1≤K≤5。
每行最后一个数字P (1≤P≤9999)表示此商品组合的优惠价。
★结果输出:将计算出的所购商品应付的最少费用输出到文件output.txt。
输入文件示例输出文件示例Input.txt offer.txt output.txt2 2 147 3 2 1 7 3 58 2 5 2 7 1 8 2 10解:设cost(a,b,c,d,e)表示购买商品组合(a,b,c,d,e)需要的最少费用。
A[k],B[k],C[k],D[k],E[k]表示第k种优惠方案的商品组合。
offer (m)是第m种优惠方案的价格。
西安交大计算方法上机报告

从而得到计算的公式:
j 1, 2,..., n 1 j 1 j l i1 i 2,3,..., n i1 11 i 1 lik ukj j i , i 1,..., n, i 2,3,.., n ij ij k 1 i 1 1 l ( lkt ti ) k i 1,..., n, i 2,3,.., n ki ki t 1 ii
计算方法 上机实习题目报告
班级:材料 s3076 班 姓名:丁明帅 学号:3113305029
1:计算 S
100000
k 1
1 ,要求误差小于 106 ,给出实现算法。 k2
【实现思路】
设当 k 值为 i 时 S-Si<10-6,则只需要
S Si
100000 1
k i
1
2
1 l32 ln 2 1 ln 3 1
4 / 46
12 22
1n 2 n
nn
因此有
1 j 2 j 0 jj 0 0
ij li1
li 2
li ,i -1
根据定义知道 L 矩阵的斜对角线上的值都为 1,且 L 矩阵的第一行与原矩阵 A 的第一行 相同,因此可以根据公式先计算 U 的第一行,然后计算 L 的第一列;以后的第 i 步先计算 U 的第 i 行, 然后计算 L 的第 i 列 (U 的第 n 行不作计算) 。 然后把最后的 U 的第 n 行计算出来。
【算法依据】
列主元高斯消元法的过程可以将方程组系数简化为系数矩阵与 b 矩阵, 从而利用方程组 对系数扩展矩阵进行消元。 在消元的过程中矩阵的行向量之间可以变换, 但列向量不能变化。 在进行压缩矩阵的求解中还需要人为的将因调整行向量所导致的列向量的变化调整回来。 Matlab 对于矩阵的处理非常容易,结合 for 循环语句可以实验本题大规模方程组的求解。
计算方法上机大作业

绪论1.用二分法求方程x3-x-1=0在[1,2]内的近似根,要求误差不超10-3. 程序:clc;cleara=1;b=2;fa=a*a*a-a-1;fb=b*b*b-b-1;c=(a+b)/2;fc=c*c*c-c-1;if fa*fb>0,break,endwhile abs(fc)>10^(-3)c=(a+b)/2;fc=c*c*c-c-1;if fb*fc>0b=c;fb=fc;elsea=c;fa=fc;endendformat longfx=fc,x=c运算结果:fx =-4.659488331526518e-05x =1.3247070312500002.证明方程f(x)=e x+10x-2=0内在[0,1]内有唯一实根,用二分法求这一实根,要求误差不超过0.5*10-2.证明:由题可知,f’(x)在[0,1]上为正,也就是f(x)单调递增,又因为f(0)<0,f(1)>0,所以在[0,1]内有唯一实根。
程序:clc;cleara=0;b=1;fa=exp(a)+10*a-2;fb=exp(b)+10*b-2;c=(a+b)/2;fc=exp(c)+10*c-2;if fa*fb>0,break,endwhile abs(fc)>0.5*10^(-2)c=(a+b)/2;fc=exp(c)+10*c-2;if fb*fc>0b=c;fb=fc;elsea=c;fa=fc;endendformat longfx=fc,x=c运算结果:fx = 0.003275341789827 x =0.090820312500000第一章11.给出概率积分dx ey xx22⎰-=π的数据表i 0 1 2 3 xi 0.46 0.47 0.48 0.49 yi0.48465550.49374520.50274980.5116683用二次插值计算,试问:(1)当x=0.472时该积分值等于多少? 程序: clear all;x=[0.46 0.47 0.48 0.49];y=[0.4846555 0.4937452 0.5027498 0.5116683]; xi=0.472;yi=interp1(x,y,xi,'linear'); yi 运算结果:yi = 0.495546120000000(2)当x 为何值时积分值等于0.5? 程序:clear all;x=[0.46 0.47 0.48 0.49];y=[0.4846555 0.4937452 0.5027498 0.5116683];yi=0.5;xi=interp1(y,x,yi,'linear');xi运算结果:xi =0.47694622748373134.构造适合下列数据表的三次插值样条S(x):x -1 0 1 2 y -1 1 3 5 y’ 6 1 程序:function [ ]=spline31(X,Y,dY,x0,m)N=size(X,2);s0=dY(1); sN=dY(2);interval=0.025;disp('x0为插值点')h=zeros(1,N-1);for i=1:N-1 h(1,i)=X(i+1)-X(i);endd(1,1)=6*((Y(1,2)-Y(1,1))/h(1,1)-s0)/h(1,1);d(N,1)=6*(sN-(Y(1,N)-Y(1,N-1))/h(1,N-1))/h(1,N-1);for i=2:N-1d(i,1)=6*((Y(1,i+1)-Y(1,i))/h(1,i)-(Y(1,i)-Y(1,i-1))/h(1,i-1))/(h(1,i)+h(1,i-1));endmu=zeros(1,N-1); md=zeros(1,N-1);md(1,N-1)=1; mu(1,1)=1;for i=1:N-2u=h(1,i+1)/(h(1,i)+h(1,i+1)); mu(1,i+1)=u;md(1,i)=1-u;endp(1,1)=2; q(1,1)=mu(1,1)/2;for i=2:N-1p(1,i)=2-md(1,i-1)*q(1,i-1); q(1,i)=mu(1,i)/p(1,i);endp(1,N)=2-md(1,N-1)*q(1,N-1);y=zeros(1,N); y(1,1)=d(1)/2;for i=2:N y(1,i)=(d(i)-md(1,i-1)*y(1,i-1))/p(1,i);endx=zeros(1,N); x(1,N)=y(1,N);for i=N-1:-1:1 x(1,i)=y(1,i)-q(1,i)*x(1,i+1);endfprintf ('M为三对角方程的解\n');M=x;fprintf ('\n');syms t;digits (m);for i=1:N-1pp(i)=M(i)*(X(i+1)-t)^3/(6*h(i))+M(i+1)*(t-X(i))^3/(6*h(i))+(Y(i)-M(i)*h(i)^2/6)*(X(i+1)-t)/h(i);+(Y(i+1)-M(i+1)*h(i)^2/6)*(x-X(i))/h(i);pp(i)=simplify(pp(i)); coeff=sym2poly(pp(i));if length(coeff)~=4tt=coeff(1:3); coeff(1:4)=0; coeff(2:4)=tt;endif x0>X(i)&&x0<X(i+1) L=i;y0=coeff(1)*x0^3+coeff(2)*x0^2+coeff(3)*x0+coeff(4); endval=X(i):interval:X(i+1);for k=1:length(val)fval(k)=coeff(1)*val(k)^3+coeff(2)*val(k)^2;+coeff(3)*val(k)+coeff(4);endif mod(i,2)==1 plot(val,fval,'r+')else plot(val,fval,'b.')endhold onclear val fvalans=sym(coeff,'d');ans=poly2sym(ans,'t');fprintf('在区间[%f,%f]内\n',X(i),X(i+1));fprintf('三次样条函数S(%d)=',i);pretty(ans);endX=[-1 0 1 3];Y=[-1 1 3 5];dY=[6 1];x0=1.2;m=5;spline31(X,Y,dY,x0,m)运行结果:在区间[-1.000000,0.000000]内三次样条函数S(1)= 3.0 t3 + 2.0 t2 + 0.66667 t + 0.66667 在区间[0.000000,1.000000]内三次样条函数S(2)= - 1.0 t3 + 2.0 t 2 - 2.3333 t + 1.0在区间[1.000000,3.000000]内三次样条函数S(3)=0.25 t 3- 1.75 t 2+ 2.5833 t + 1.9167第二章PPT题目:编写一通用型复化辛甫生公式。
奥鹏西安交通大学课程考试《运筹学》参考资料答案.doc

奥鹏西安交通⼤学课程考试《运筹学》参考资料答案.doc 西安交通⼤学课程考试复习资料单选题1.从甲市到⼄市之间有-公路⽹络,为了尽快从甲市驱车赶到⼄市,应借⽤()A.树的逐步⽣成法B.求最⼩技校树法C.求最短路线法D.求最⼤流量法答案: C2.⼯序A是⼯序B的紧后⼯序,则错误的结论是A.⼯序B完⼯后⼯序A才能开⼯B.⼯序A完⼯后⼯序B才能开⼯C.⼯序B是⼯序A的紧前⼯序D.⼯序A是⼯序B的后续⼯序答案: B3.线性规划的求解中,⽤最⼩⽐值原则确定换出变量,⽬的是保持解的可⾏性。
()A.正确B.错误C.不⼀定D.⽆法判断答案: A4.⽤图解法求解⼀个关于最⼤利润的线性规划问题时,若其等利润线与可⾏解区域相交,但不存在可⾏解区域最边缘的等利润线,则该线性规划问题( )。
A.有⽆穷多个最优解B.有可⾏解但⽆最优解C.有可⾏解且有最优解D.⽆可⾏解答案: B5.线性规划的图解法中,⽬标函数值的递增⽅向与()有关?A.约束条件B.可⾏域的范围C.决策变量的⾮负性D.价值系数的正负最好挑选( )为调整格。
A.WB格B.WC格C.YA格D.XC格答案: A7.线性规划的图解法中,⽬标函数值的递增⽅向与()有关?A.约束条件B.可⾏域的范围C.决策变量的⾮负性D.价值系数的正负答案: D8.⽤运筹学解决问题时,要对问题进⾏()A.分析与考察B.分析和定义C.分析和判断D.分析和实验答案: B9.影⼦价格的经济解释是()A.判断⽬标函数是否取得最优解B.价格确定的经济性C.约束条件所付出的代价D.产品的产量是否合理答案: C10.求解线性规划模型时,引⼊⼈⼯变量是为了()A.使该模型存在可⾏解B.确定⼀个初始的基可⾏解C.使该模型标准化D.其他均不正确答案: B11.⼀般讲,对于某⼀问题的线性规划与该问题的整数规划可⾏域的关系存在()A.前者⼤于后者B.后者⼤于前者12.影⼦价格的经济解释是()A.判断⽬标函数是否取得最优解B.价格确定的经济性C.约束条件所付出的代价D.产品的产量是否合理答案: C13.在⼀个运输⽅案中,从任⼀数字格开始,( )⼀条闭合回路。
西北工业大学智慧树知到“计算机科学与技术”《计算方法》网课测试题答案卷3

西北工业大学智慧树知到“计算机科学与技术”《计算方
法》网课测试题答案
(图片大小可自由调整)
第1卷
一.综合考核(共10题)
1.利用待定系数法可以得出各种求积公式,而且可以具有尽可能高的代数精度。
()
A.正确
B.错误
2.判断参数值是否正确{图}。
()
A.正确
B.错误
3.牛顿迭代法的基本思想是将非线方程f(x)=0逐步转化为线性议程来求解。
()
A.正确
B.错误
4.雅可比方法的主要特点是什么()
A.精度高
B.算法稳定
C.稀疏性
D.求得的特征向量正交性好
5.议程的近似方法有()
A.迭代法
B.牛顿法
C.弦截法
D.二分法
6.{图}1
A.正确
B.错误
7.直接法是在理论上没有舍入误差的前提下经过有限步运算即可得到方程组的精确解。
() A.正确
B.错误
8.{图}1
A.D
B.C
C.B
D.A
9.列主元素消元法不是直接法中常用的方法。
()
A.正确
B.错误
10.乘幂法主要是用来求矩阵的主特征值(按模最大的特征值)及相应的特征向量。
()
A.正确
B.错误
第1卷参考答案
一.综合考核
1.参考答案:A
2.参考答案:A
3.参考答案:A
4.参考答案:ABD
5.参考答案:ABCD
6.参考答案:A
7.参考答案:A
8.参考答案:D
9.参考答案:B
10.参考答案:A。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
计算方法A 上机大作业1. 共轭梯度法求解线性方程组算法原理:由定理3.4.1可知系数矩阵A 是对称正定矩阵的线性方程组Ax=b 的解与求解二次函数1()2TT f x x Ax b x =- 极小点具有等价性,所以可以利用共轭梯度法求解1()2TT f x x Ax b x =-的极小点来达到求解Ax=b 的目的。
共轭梯度法在形式上具有迭代法的特征,在给定初始值情况下,根据迭代公式:(1)()()k k k k x x d α+=+产生的迭代序列(1)(2)(3)x x x ,,,... 在无舍入误差假定下,最多经过n 次迭代,就可求得()f x 的最小值,也就是方程Ax=b 的解。
首先导出最佳步长k α的计算式。
假设迭代点()k x 和搜索方向()k d 已经给定,便可以通过()()()()k k f x d φαα=+的极小化()()min ()()k k f x d φαα=+来求得,根据多元复合函数的求导法则得:()()()'()()k k T k f x d d φαα=∇+令'()0φα=,得到:()()()()k T k k k T k r d d Adα= ,其中()()k k r b Ax =-然后确定搜索方向()k d 。
给定初始向量(0)x 后,由于负梯度方向是函数下降最快的方向,故第一次迭代取搜索方向(0)(0)(0)(0)()dr f x b Ax ==-∇=- 。
令(1)(0)00x x d α=+其中(0)(0)0(0)(0)T T r d d Adα=。
第二次迭代时,从(1)x 出发的搜索方向不再取(1)r ,而是选取(1)(1)(0)0dr d β=+,使得(1)d 与(0)d 是关于矩阵A 的共轭向量,由此可求得参数0β:(1)(0)0(0)(0)T T r Ad d Adβ=-然后从(1)x 出发,沿(1)d 进行搜索得到(2)(1)(1)1x x d α=+设已经求出(1)()()k k k k x x d α+=+,计算(1)(1)k k r b Ax ++=-。
令(1)(1)()k k k k dr d β++=+,选取k β,使得(1)k d +和()k d 是关于A 的共轭向量,可得:(1)()()()k T k k k T k r Ad d Adβ+=-具体编程计算过程如下:(1) 给定初始近似向量(0)x 以及精度0ε> ; (2) 计算(0)(0)r b Ax =-,取(0)(0)d r =; (3) For k=0 to n-1 do(i )()()0()()k T k k T k r d d Adα=;(ii )(1)()()k k k k xx d α+=+;(iii )(1)(1)k k r b Ax ++=-;(iv )若(1)k r ε+≤或k=n-1,则输出近似解(1)k x + ,停止;否则,转(v );(v )2(1)22()2k k k r rβ+= ;(vi )(1)(1)()k k k k d r d β++=+;End do程序框图:程序使用说明:本共轭梯度法求解线性方程的程序直接打开matlab运行,在求解线性方程组Ax=b(A是对称正定矩阵)的时候,直接运行程序Gongetidufa,输入A,b的值,虽然该函数是通用的,但是对于矩阵A和向量b的输入,需要使用者根据A和b的特点自行输入。
算例3.4.2计算结果:对99页例题3.4.2,运行程序Gongetidufa将矩阵A,b读入系统程序如下:clear allclcA=input('请输入A的值'); %输入n阶正定矩阵A的值b=input('请输入b的值:'); %输入b的值n=size(A,1); %求出矩阵A的行数x=zeros(n,1); %给定x的初始值e=10^(-12); %给出精度r=b-A*x;d=r;for(i=1:n)a=norm(r,2)^2/(d'*A*d); %求最佳步长 x=x+a*d; j=r;r=b-A*x;if(norm(r)<=e||i==n) break; elseB=norm(r,2)^2/norm(j,2)^2; d=r+B*d; end endx %输出最终的x 的结果 计算结果:x=[1;1;1]2.三次样条差值算法原理(三次样条插值函数的导出): (i).导出在子区间[]1,i i x x - 上的S(x)的表达式由于S(x)的二阶导数连续,设S(x)再节点i x 处的二阶导数值S ’’(xi)=Mi ,其中Mi 为未知的待定参数。
由S(x)是分段的三次多项式知,S ’’(x)是分段线性函数,S ’’(x)在子区间[]1,i i x x -上可表示为1111111''(),i i i ii i i i i i i i i ii ix x x x S x M M x x x x x x x x M M x x x h h ---------=+----=+≤≤其中hi=xi-x(i-1),对上式两次积分得到()()331111()66()(),i i i i iii i i i i ix x x x S x M M h h b x x c x x x x x ------=++-+-≤≤由插值条件11(),()i i i i S x y S x y --== 得到2211()/,()/66i i i i i i i i i i h h b y M h c y M h --=-=-将i i b c 和 代入()S x 可得()()3321111211()()/()666()/(),6ii i i i i i i i i iii i i i i i x x x x h S x M M y M h x x h h hy M h x x x x x --------=++--+--≤≤(ii).建立参数i M 的方程组 对S(x)求导可得()()2211111'()()/22,6ii i i i i iiii i i i ix x x x S x M M y y h h h M M h x x x -------=-++---≤≤上式中令i x x = 得S(x)在xi 处的左导数'()i S x - ,令1i x x -=得到右导数'()i S x +,因为S(x)在内节点xi 处一阶导数连续,所以''()()i i S x S x -+=,进一步推导可得112,1,2,3,...,1i i i i i i M M M d i n μλ-+++==-其中1ii i i h h h μ+=+,111i i i i i h h h λμ++==-+,1111116()6[,,],1,2,...,1i i i i i i i i i i i iy y y y d f x x x i n h h h h +--+++--=-==-+上式为三弯矩方程组,因为三弯矩方程组只有n-1个方程,不能确定n+1个未知量Mi ,所以需要再增加两个方程,由边界条件确定。
第一种边界条件:此时已知''()''()f a f b 和 .不妨取0''()M f a =,''()n M f b =,这时三弯矩方程组化为:1121101112111222,3,...,22i i i i i in i n n n n M M d M M M M d i n M M d Mλμμλμλ-+-----+=-⎧⎪++==-⎨⎪+=-⎩ 以上方程组系数矩阵式严格三对角占优矩阵,可用追赶法求解。
求出(1,2,...,1)i M i n =- 后,代入S(x)可得三次样条插值函数的数学表达式。
第二种边界条件:已知'()'()f a f b 和。
记0''()''()n y f a y f b ==,,则有00'()''()'n n S x y S x y +-==,所以:1011101011',',3663n n n n n n n ny y h h y y h hM M y M M y h h ------+=++= 即0102M M d +=12n n n M M d -+=其中10001116(')6(')n n n n n ny y d y h h y y d y h h --=--=-所以得到第二种边界条件下的三弯矩方程组:0101112,2,1,2,3,...,1,2i i i i i i n n n M M d M M M d i n M M dμλ-+-+=⎧⎪++==-⎨⎪+=⎩ 该方程组系数矩阵是严格三对角占优矩阵,可用追赶法求解,具体追赶法的求解过程见《数值分析》教材。
第三种边界条件:周期型边界条件.已知()y f x = 是以0n T b a x x =-=- 为周期的周期函数,则由周期性可知,01101111,,,,n n n n n y y y y M M M M h h +++===== ,这时,可以将点n x 看成内点,则方程组对i=n 也成立,既有112n n n n n n M M M d μλ-+++= ,也即112n n n n n M M M d λμ-++=,其中11116()n n n n n ny y y y d h h h h ---=-+ 于是三弯矩方程组化为1121111112,2,2,3,4,....,1,2.n i i i i i i n n n n n M M M d M M M d i n M M M d λμμλλμ-+-++=⎧⎪++==-⎨⎪++=⎩ 该方程组可用matlab 直接解出。
程序框图如下:程序使用说明:本程序是求解137页例题4.6.1的运行结果,通过程序便可求得M ,然后根据3321111211()()()()666(),6i i i ii i i i i i iii i i i i ix x x x h x x S x M M y M h h h h x x y M x x x h ---------=++--+-≤≤便可得到,在1i i x x x -≤≤ 上的三次样条插值函数()S x ,进而得到整个区间上的三次样条差值函数()S x 。
算例计算结果:137页例题4.6.1的计算实习1、打开matlab 运行Sanciyangtiao 程序2、自行输入x 和y 的节点值3、得出计算结果3.龙贝格积分法对于复化梯形求积公式,取积分近似值2221()()41n n n n I f T T T T ≈+-=- 对复化辛普森求积公式4(4)()(),2880n b a I f S h f a b ηη--≈-≤≤4(4)211()()(),28802n b a h I f S f a b ηη--≈-≤≤ 因为(4)(4)1()()ff ηη=,所以上述两式相除得2()16()nnI f S I f S -≈-所以,22221()()41n n n n I f S S S S ≈+-=-同理,22231()()41n n n n I f C C C C ≈+-=-对2n T ,2n S 和2n C 分析可得222222231()411()411()41n n n n n n n n n n n n S T T T C S S S R C C C ⎧=+-⎪-⎪⎪=+-⎨-⎪⎪=+-⎪-⎩龙贝格积分算法如下:1111111121122122222222232222[()()],21((21)),0,1,2...,2221(),411(),0,1,2...,411(),41kk k k k k k k k k k k k k k k k i b a T f a f b b a b aT T f a i k S T T T C S S S k R C C C +++++++++=-=+--=++-=⎧=+-⎪-⎪⎪=+-=⎨-⎪⎪=+-⎪-⎩∑如果122,k k R R ε+-< 则取12[]k I f R +≈,否则,继续计算直到满足条件。