上海大学课程实验报告-数值代数

合集下载

数值代数实验报告

数值代数实验报告

数值代数实验报告数值代数实验报告引言:数值代数是一门研究数值计算方法和算法的学科,它在科学计算和工程应用中起着重要的作用。

本实验报告旨在通过实际的数值计算问题,探讨数值代数的应用和效果。

实验一:线性方程组求解线性方程组求解是数值代数中的一个重要问题。

在实验中,我们使用了高斯消元法和LU分解法两种求解线性方程组的方法,并对比了它们的效果。

首先,我们考虑一个3×3的线性方程组:2x + 3y - z = 54x - 2y + 2z = 1x + y + z = 3通过高斯消元法,我们将该方程组转化为上三角形式,并得到解x=1, y=2, z=0。

而通过LU分解法,我们将该方程组分解为LU两个矩阵的乘积,并得到相同的解。

接下来,我们考虑一个更大的线性方程组,例如10×10的方程组。

通过比较高斯消元法和LU分解法的运行时间,我们可以发现LU分解法在处理大规模方程组时更加高效。

实验二:特征值与特征向量计算特征值与特征向量计算是数值代数中的另一个重要问题。

在实验中,我们使用了幂法和QR方法两种求解特征值与特征向量的方法,并对比了它们的效果。

首先,我们考虑一个3×3的矩阵:1 2 34 5 67 8 9通过幂法,我们可以得到该矩阵的最大特征值为15.372,对应的特征向量为[0.384, 0.707, 0.577]。

而通过QR方法,我们也可以得到相同的结果。

接下来,我们考虑一个更大的矩阵,例如10×10的矩阵。

通过比较幂法和QR 方法的运行时间,我们可以发现QR方法在处理大规模矩阵时更加高效。

实验三:奇异值分解奇异值分解是数值代数中的一种重要技术,它可以将一个矩阵分解为三个矩阵的乘积,从而实现数据降维和信息提取的目的。

在实验中,我们使用了奇异值分解方法,并通过实际的数据集进行了验证。

我们选取了一个包含1000个样本和20个特征的数据集,通过奇异值分解,我们将该数据集分解为三个矩阵U、S和V的乘积。

数值代数上机实验报告

数值代数上机实验报告

数值代数上机实验报告试验项目名称:平方根法与改进平方根法实验内容:先用你熟悉的计算机语言将平方根法和改进平方根法编写成通用的子程序,然后用你编写的程序求解对称正定方程组Ax=b,其中,A=[101 10 1…1 10 11 10]100*100b随机生成,比较计算结果,评论方法优劣。

实验要求:平方根法与改进的平方根的解法步骤;存储单元,变量名称说明;系数矩阵与右端项的生成;结果分析。

实验报告姓名:罗胜利班级:信息与计算科学0802 学号:u200810087 实验一、平方根法与改进平方根法先用你所熟悉的计算机语言将平方根法和改进的平方根法编成通用的子程序,然后用你编写的程序求解对称正定方程组AX=b,其中系数矩阵为40阶Hilbert矩阵,即系数矩阵A的第i行第j列元素为=,向量b的第i个分量为=.平方根法函数程序如下:function [x,b]=pingfanggenfa(A,b)n=size(A);n=n(1);x=A^-1*b; %矩阵求解disp('Matlab自带解即为x');for k=1:nA(k,k)=sqrt(A(k,k));A(k+1:n,k)=A(k+1:n,k)/A(k,k);for j=k+1:n;A(j:n,j)=A(j:n,j)-A(j:n,k)*A(j,k);endend %Cholesky分解for j=1:n-1b(j)=b(j)/A(j,j);b(j+1:n)=b(j+1:n)-b(j)*A(j+1:n,j);endb(n)=b(n)/A(n,n); %前代法A=A';for j=n:-1:2b(j)=b(j)/A(j,j);b(1:j-1)=b(1:j-1)-b(j)*A(1:j-1,j);endb(1)=b(1)/A(1,1); %回代法disp('平方根法的解即为b');endfunction [x]=ave(A,b,n) %用改进平方根法求解Ax=b L=zeros(n,n); %L为n*n矩阵D=diag(n,0); %D为n*n的主对角矩阵S=L*D;for i=1:n %L的主对角元素均为1L(i,i)=1;for i=1:nfor j=1:n %验证A是否为对称正定矩阵if (eig(A)<=0)|(A(i,j)~=A(j,i)) %A的特征值小于0或A非对称时,输出wrong disp('wrong');break;endendendD(1,1)=A(1,1); %将A分解使得A=LDL Tfor i=2:nfor j=1:i-1S(i,j)=A(i,j)-sum(S(i,1:j-1)*L(j,1:j-1)');L(i,1:i-1)=S(i,1:i-1)/D(1:i-1,1:i-1);endD(i,i)=A(i,i)-sum(S(i,1:i-1)*L(i,1:i-1)');endy=zeros(n,1); % x,y为n*1阶矩阵x=zeros(n,1);for i=1:ny(i)=(b(i)-sum(L(i,1:i-1)*D(1:i-1,1:i-1)*y(1:i-1)))/D(i,i); %通过LDy=b解得y的值endfor i=n:-1:1x(i)=y(i)-sum(L(i+1:n,i)'*x(i+1:n)); %通过L T x=y解得x的值改进平方根法函数程序如下:function b=gaijinpinfanggenfa(A,b)n=size(A);n=n(1);v=zeros(n,1);for j=1:nfor i=1:j-1v(i)=A(j,i)*A(i,i);endA(j,j)=A(j,j)-A(j,1:j-1)*v(1:j-1);A(j+1:n,j)=(A(j+1:n,j)-A(j+1:n,1:j-1)*v(1:j-1))/A(j,j);end %LDL'分解B=diag(A);D=zeros(n);for i=1:nD(i,i)=B(i);A(i,i)=1;EndA=tril(A); %得到L和Dfor j=1:n-1b(j)=b(j)/A(j,j);b(j+1:n)=b(j+1:n)-b(j)*A(j+1:n,j);endb(n)=b(n)/A(n,n); %前代法A=D*(A');for j=n:-1:2b(j)=b(j)/A(j,j);b(1:j-1)=b(1:j-1)-b(j)*A(1:j-1,j);endb(1)=b(1)/A(1,1); %回代法disp('改进平方根法解得的解即为b');end调用函数解题:clear;clc;n=input('请输入矩阵维数:');b=zeros(n,1);A=zeros(n);for i=1:nfor j=1:nA(i,j)=1/(i+j-1);b(i)=b(i)+1/(i+j-1);endend %生成hilbert矩阵[x,b]=pingfanggenfa(A,b) b=gaijinpinfanggenfa(A,b)运行结果:请输入矩阵维数:40Warning: Matrix is close to singular or badly scaled. Results may be inaccurate. RCOND = 6.570692e-020. > In pingfanggenfa at 4In qiujie at 10Matlab自带解即为x平方根法的解即为bx =1.60358.96850.85621.01950.9375-50.2500-3.0000-16.000024.0000-49.5000-30.000039.000022.0000-64.0000-12.00002.000010.2500-10.5000-1.0000-10.875083.000046.0000-98.0000-69.000068.000021.0000-50.7188-8.7500-8.0000 112.0000 6.0000 -68.7500 22.000044.0000 -28.0000 8.0000 -44.000012.0000b =1.0e+007 *0.0000-0.00000.0001-0.0004-0.00140.0424-0.29801.1419-2.73354.2539-4.30182.7733-1.19890.5406-0.36880.32850.4621-0.25130.05650.0000-0.00510.0071-0.0027-0.0031-0.00190.00090.0002-0.0002-0.00060.00040.0001-0.00020.00010.0000-0.00000.0000-0.0000-0.0000改进平方根法解得的解即为bb =1.0e+024 *0.0000-0.00000.0001-0.0012-0.0954 0.4208 -1.2101 2.0624 -1.0394 -3.3343 6.2567 -0.2463 -7.45942.80303.6990 0.7277 -1.7484 -0.4854 -3.6010 0.2532 5.1862 1.4410 0.8738 -4.5654 1.0422 4.0920 -2.7764 -2.2148 -0.8953 0.3665 4.8967 1.0416 0.1281-1.1902-2.83348.4610-3.6008实验二、利用QR分解解线性方程组:利用QR分解解线性方程组Ax=b,其中A=[16 4 8 4;4 10 8 4;8 8 12 10;4 4 10 12];b=[32 26 38 30];求解程序如下:定义house函数:function [v,B]=house(x)n=length(x);y=norm(x,inf);x=x/y;Q=x(2:n)'*x(2:n);v(1)=1;v(2:n)=x(2:n);if n==1B=0;elsea=sqrt(x(1)^2+Q);if x(1)<=0v(1)=x(1)-a;elsev(1)=-Q/(x(1)+a);endB=2*v(1)^2/(Q+v(1)^2);endend进行QR分解:clear;clc;A=[16 4 8 4;4 10 8 4;8 8 12 10;4 4 10 12]; b=[32 26 38 30];b=b';x=size(A);m=x(1);n=x(2);d=zeros(n,1);for j=1:n[v,B]=house(A(j:m,j));A(j:m,j:n)=(eye(m-j+1)-B*(v')*v)*A(j:m,j:n); d(j)=B;if j<m< p="">A(j+1:m,j)=v(2:m-j+1);endend %QR分解R=triu(A); %得到R D=A;I=eye(m,n);Q=I;for i=1:nD(i,i)=1;endH=tril(D);M=H';for i=1:nN=I-d(i)*H(1:m,i)*M(i,1:m);Q=Q*N;end %得到Qb=(Q')*b; %Q是正交阵for j=n:-1:2b(j)=b(j)/R(j,j);b(1:j-1)=b(1:j-1)-b(j)*R(1:j-1,j);endb(1)=b(1)/R(1,1); %回带法运行结果如下:R =18.7617 9.8072 15.7769 11.08640 9.9909 9.3358 7.53410 0 5.9945 9.80130 0 0 -0.5126Q =0.8528 -0.4368 -0.2297 -0.17090.2132 0.7916 -0.4594 -0.34170.4264 0.3822 0.2844 0.76890.2132 0.1911 0.8095 -0.5126b=1.000000000000001.000000000000010.9999999999999881.00000000000001实验三、Newton下山法解非线性方程组:3x-cos(yz)-=0,-81+sinz+1.06=0,exp(-xy)+20z+=0;要求满足数值解=满足或.定义所求方程组的函数:Newtonfun.mfunction F = Newtonfun(X)F(1,1)=3*X(1)-cos(X(2)*X(3))-1/2;F(2,1)=X(1)^2-81*(X(2)+0.1)^2+sin(X(3))+1.06;F(3,1)=exp(-X(1)*X(2))+20*X(3)+(10*pi-3)/3;End向量求导:Xiangliangqiudao.mfunction J=xiangliangqiudao()syms x y zX=[x,y,z];F=[3*X(1)-cos(X(2)*X(3))-1/2;X(1)^2-81*(X(2)+0.1)^2+sin(X(3))+1.06;exp(-X(1)*X(2))+20*X(3)+(10*pi-3)/3];J=jacobian(F,[x y z]);End代值函数:Jacobi.mfunction F=Jacobi(x)F=[ 3,x(3)*sin(x(2)*x(3)), x(2)*sin(x(2)*x(3));2*x(1), -162*x(2)-81/5,cos(x(3));-x(2)/exp(x(1)*x(2)),-x(1)/exp(x(1)*x(2)),20];End方程组求解:format long; %数据表示为双精度型X1=[0,0,0]';eps=10^(-8);k=1;i=1;X2=X1-Jacobi(X1)^(-1)*Newtonfun(X1);while (norm(subs(X2-X1,pi,3.1415926),2)>=eps)&&(norm(Newtonfun(X1),2)>=eps) if norm(Newtonfun(X2),2)<="" p="">X1=X2;B=inv(Jacobi(X2));C=Newtonfun(X2);X2=X2-B*C;i=i+1;elsev=1/(2^k); %引入下山因子X1=X2;B=inv(Jacobi(X2));C=Newtonfun(X2);X2=X2-v*B*C;k=k+1;endendj=i+k-1 %迭代次数X=X2 %输出结果运行结果如下:j =5X =0.500000000000000 -0.000000000000000 -0.523598775598299</m<>。

实验报告——常微分方程的数值解法

实验报告——常微分方程的数值解法

实验报告实验项目名称常微分方程的数值解法实验室数学实验室所属课程名称微分方程数值解实验类型上机实验实验日期2013年3月11日班级10信息与计算科学学号2010119421姓名叶达伟成绩实验概述:【实验目的及要求】运用不同的数值解法来求解具体问题,并通过具体实例来分析比较各种常微分方程的数值解法的精度,为以后求解一般的常微分方程起到借鉴意义。

【实验原理】各种常微分方程的数值解法的原理,包括Euler法,改进Euler法,梯形法,Runge-Kutta方法,线性多步方法等。

【实验环境】(使用的软硬件)Matlab软件实验内容:【实验方案设计】我们分别运用Euler法,改进Euler法,RK方法和Adams隐式方法对同一问题进行求解,将数值解和解析解画在同一图像中,比较数值解的精度大小,得出结论。

【实验过程】(实验步骤、记录、数据、分析)我们首先来回顾一下原题:对于给定初值问题:1. 求出其解析解并用Matlab画出其图形;2. 采用Euler法取步长为0.5和0.25数值求解(2.16),并将结果画在同一幅图中,比较两者精度;3. 采用改进Euler法求解(2.16),步长取为0.5;4. 采用四级Runge-Kutta法求解(2.16),步长取为0.5;5. 采用Adams四阶隐格式计算(2.16),初值可由四级Runge-Kutta格式确定。

下面,我们分五个步骤来完成这个问题:步骤一,求出(2.16)式的解析解并用Matlab 画出其图形; ,用Matlab 做出函数在上的图像,见下图:00.51 1.52 2.53 3.54 4.550.511.522.533.5x 1015y=exp(1/3 t 3-1.2t)exact solution图一 初值问题的解析解的图像步骤二,采用Euler 法取步长为0.5和0.25数值求解(2.16),并将结果画在同一幅图中,比较两者精度;我们采用Euler 法取步长为0.5和0.25数值求解,并且将数值解与解析解在一个图中呈现,见下图:00.51 1.52 2.53 3.54 4.550.511.522.533.5x 1015Numerical solution of Euler and exact solutionexact solution h=0.5h=0.25图二 Euler 方法的计算结果与解析解的比较从图像中不难看出,采用Euler 法取步长为0.5和0.25数值求解的误差不尽相同,也就是两种方法的计算精度不同,不妨将两者的绝对误差作图,可以使两种方法的精度更加直观化,见下图:00.51 1.52 2.53 3.54 4.550.511.522.533.5x 1015Absolute error of numerical solution and exact solutionh=0.5h=0.25图三 不同步长的Euler 法的计算结果与解析解的绝对误差的比较 从图像中我们不难看出,步长为0.25的Euler 法比步长为0.5的Euler 法的精度更高。

《数值计算》实验报告

《数值计算》实验报告

《数值计算》实验报告第一部分:简答题(请简要回答以下问题,每小题字数不少于200字)1、Matlab变量命名有什么要求?以下变量名是否合法?对不合法的变量名说明理由。

abcd-2xyz_33chan NaN ABCDefgh2、插值、拟合、回归这三种方法是用来解决什么问题的?面对一组数据,如何选择用什么方法?3、数值积分的主要思想是什么?常用的数值积分公式有哪几个?4、请结合自己的学习,举例说明《数值计算》课程中所学方法在解决实际问题中是如何应用的。

第二部分:基础题(请完成以下问题,要求给出程序语句及计算结果,用截图方式附在各题目下方)1、已知点(1,3.0),(2,3.7),(5,3.9),(6,4.2),(7,5.7),(8,6.6),(10,7.1),(13,6.7), (17,4.5),绘出经过这些点的函数曲线图形,并给出曲线方程。

答:采取三次样条插值法,九个输入数据分成八段,每一段就是一个三次函数。

这八段的函数形式为y = a0 + a1*x + a2*x^2 + a3*x^3,每个分段函数的参数构成下图所示的coefs 矩阵。

2、在我国某海域测得海洋不同深度处的水温如表1所示,求水深为800m和1500m处的温度。

答:采取线性插值法求得800m和1500m处的温度3、求解方程组⎪⎪⎩⎪⎪⎨⎧=-++=--=-++=++56533332821w z y x w y x w z y x z y x ,请至少使用两种方法求解,并对这两种方法的计算结果进行说明。

高斯消元法LU 分解QR分解Jacobi迭代法使用 Jacobi 迭代法无法求出结果,表示迭代的过程中不收敛4、计算积分dx eI x ⎰-=1022,精度为10-6。

被积函数总共调用 13 次,求得积分值为 0.85565、求方程t et t f t5.0)(sin )(1.02-⋅=-在[0.5,1]内的根。

6、求解微分方程0)1(22=+'--''y y y y ,0)0(,1)0(,300='=≤≤y y x ,绘出解函数的图形。

数值代数(第13周)实验报告()

数值代数(第13周)实验报告()
[x,k]=fastest(A,b,eps);
运行结果:
x =
-4.0000
3.0000
2.0000
k =
17
(2)第二题的第一个方程组
A=[1 0.4 0.4;0.4 1 0.8;0.4 0.8 1];
b=[1 2 3]';
eps=0.00001;
[x,k]=fastest(A,b,eps);
运行结果:
(本次实验进行总结,包括但不限于对算法的理解、编程实验中遇到的问题、实验中所取得的经验等;如程序未能通过,应分析错误原因.)
通过不断试验最速下降法和范数norm(x-x0)取值范围有关.前两个方程组当eps取得越小误差越小。上面的运算结果则是最接近值的迭代次数
x =
-0.1351
-1.0811
3.9189
k =
51
(3)第二题的第二个方程组
A=[1 2 -2;1 1 1;2 2 1];
b=[1 1 1]';
eps=0.1;
[x,k]=fastest(A,b,eps);
运行结果:x =
-3.0000
3.0000
1.0000
k =
5
4实验总结【必填栏目】
得分:
x0=x;
ifk>=m
disp('迭代次数太多,可能 Nhomakorabea收敛!');
return;
end
end
x
k
3.2实验(或测试)代码和结果
得分:
用最速下降法求解课本209页习题1、2中的三个方程组.
(1)第一题的方程组
A=[5 2 1;-1 4 2;2 -3 10];

上海大学数值方法报告

上海大学数值方法报告

2013-2014学年冬季学期数值方法实验报告组别第X组学号1212XXX姓名XXXX指导老师XXXX完成日期201X.XXX实验一一、题目P311.根据习题12和习题13构造算法和MATLAB 程序,以便精确计算所有情况下的二次方程的根,包括ac b b 42-≈的情况。

2.参照例1.25,对下列3个序列求序列∞=⎭⎬⎫⎩⎨⎧121n n ,请计算出前10个数值近似值。

构造类似表1.4、表1.5以及图1.8至图1.10的输出。

(a) 994.00=r ;121-=n n r r ,初始误差为0.00 2其中n=1,2,…(b) 10=q ,497.01=q ,2125---=n n n q q q ,初始误差为0.003,其中n=2,3,…二、代码第一题:第二题:三、结果1、2、四、总结本次作业的目的在于熟悉我们对Matlab基本的操作。

第一题是对if…else…的应用。

第二题是画图和格式输出。

实验二一、题目P401. 使用程序2.1求解下面每个函数的不动点(尽可能多)近似值,答案精确到小数点后12位。

同时,构造每个函数和直线y=x 来显示所有不动点。

(a ) 223)(235+--=x x x x g(b ) ))cos(sin()(x x g = (c ) )15.0()(2+-=x in x x g (d ))cos()(x x x x g -=P493. 修改程序2.2和程序2.3,使得输出分别类似于表2.1和表2.2的矩阵(即矩阵的第一行应当为[0 0a 0c 0b )(0c f ] )二、代码第一题:此题包含了文件fixpt.m 、plotfixpt.m 、sqrtm.m 、main1.m fixpt.mplotfixpt.msqrtm.mmain1.mfunction mainfigure(1)hold offgrid onaxis([-2 2.5 -2 2.5]);axis squarehold ong=inline('x.^5-3*x.^3-2*x.^2+2');x=[-2:0.01:2.5];plot(x,g(x),'r');plot(x,x,'g');f1=inline('sqrtm(3*x.^3+2*x.^2+x-2,5)');plot(x,f1(x),'b');f2=inline('sqrtm((x.^5-2*x.^2-x+2)/3,3)');plot(x,f2(x),'c');f3=inline('sqrt((x.^5-3*x.^3-x+2)./2)');plot(x,f3(x),'m');text(-0.9,-0.9,'y=x');text(-0.4,2.1,'g(x)=x^5-3x^3-2*x^2+2');text(-1.3,-1.6,'f1(x)=(3x^3+2x^2+x-2)^0^/^5');text(-0.9,0.3,'f2(x)=((x^5-2x^2-x+2)/3)^1^/^3');text(-1,1.6,'f3(x)=((x^5-3x^3-x+2)/2)^1^/^2');%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% figure(2)hold offgrid onaxis([-2 2.5 -2 2.5]);axis squarehold on[k,p,err,P]=plotfixpt(f1,0.9,0.00001,100,'k');[k,p,err,P]=plotfixpt(f1,0.6,0.00001,100,'k');plot(x,g(x),'r');plot(x,x,'g');plot(x,f1(x),'b');text(-0.9,-0.9,'y=x');text(-0.4,2.1,'g(x)=x^5-3x^3-2*x^2+2');text(-1.3,-1.6,'f1(x)=(3x^3+2x^2+x-2)^0^/^5');%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% figure(3)hold offgrid onaxis([-2 2.5 -2 2.5]);axis squarehold on[k,p,err,P]=plotfixpt(f2,-1,0.00001,100,'k');plot(x,g(x),'r');plot(x,x,'g');plot(x,f2(x),'c');text(-0.9,-0.9,'y=x');text(-0.4,2.1,'g(x)=x^5-3x^3-2*x^2+2');text(-0.9,0.3,'f2(x)=((x^5-2x^2-x+2)/3)^1^/^3');main2function mainfigure(1)hold offgrid onaxis([0 2 0 2]);axis squarehold ong=inline('cos(sin(x))');x=[0:0.01:2];plot(x,g(x),'r');plot(x,x,'g');[k,p,err,P]=plotfixpt(g,1.5,0.00001,100,'k');第二题function maing=inline('x.*sin(x)-1');[c,err,a,b,yc]=bisect(g,0,2,0.0000001);fprintf(' x sin(x)-1=0\n');fprintf('--------------------------------------------------------------------------------\n');fprintf(' k | ak | ck | bk | f(ck) \n');fprintf('---|------------------|------------------|-------------------|-------------------\n');maxk=length(yc);for k=1:maxkfprintf('%2d | %13.8f | %13.8f | %13.8f | %13.8f\n',k-1,a(k),c(k),b(k),yc(k));endfprintf('---|------------------|------------------|-------------------|-------------------\n');%g=inline('x.*sin(x)-1');[c,err,a,b,yc]=regula(g,0,2,0.0000001,0.0000001,30);fprintf(' x sin(x)-1=0\n');fprintf('--------------------------------------------------------------------------------\n');fprintf(' k | ak | ck | bk | f(ck) \n');fprintf('---|------------------|------------------|-------------------|-------------------\n'); maxk=length(yc);for k=1:maxkfprintf('%2d | %13.8f | %13.8f | %13.8f | %13.8f\n',k-1,a(k),c(k),b(k),yc(k));endfprintf('---|------------------|------------------|-------------------|-------------------\n');三、结果第一题:第二题:四、总结本次实验是解方程。

数值代数qr方法实验报告

数值代数qr方法实验报告

一、引言很多情况下在工程中抽象出来的数学方程组是超定的,没有精确解,这样就需要找一个在某种意义下最接近精确解的解。

设A 是m*n 的实矩阵,22A (A )T x bQ x b -=-,这样2min A y b -就等价与2min (A )T Q x b -,为方便求解,需要A T Q 是上三角矩阵,这样引入QR 分解就比较求解方便。

二、豪斯霍尔德变换豪斯霍尔德变换(Householder transformation )又称初等反射(Elementary reflection ),最初由A.C Aitken 在1932年提出,Alston Scott Householder 在1958年指出了这一变换在数值线性代数上的意义。

这一变换将一个向量变换为由一个超平面反射的镜像,是一种线性变换。

其变换矩阵被称作豪斯霍尔德矩阵,在一般内积空间中的类比被称作豪斯霍尔德算子,超平面的法向量被称作豪斯霍尔德向量。

三、理论依据householder 变换:任意的一个向量x ,经过一个正交变换2**T H I W W =- 后总可以变为一个与之范数相等的另一个向量Hx 。

如上图中所示,记v Hx x =-,2/w v v =,2**TH I W W =-上述H 既为所要求的householder 变换。

具体操作时将需要变化的矩阵的每一列当做一个向量,第一列变为除了第一个元素都为0的向量…第i 列变为除了前i 个元素都为0的向量...为了保证在每次变化时不改变已经变好的0元素,第i 次只变化每列第i 个元素到第n 个元素,每个变化矩阵的形式是[I 0;0 H]。

算法如下:qr_Houserhoilder.m function [Q,R]=qr_Houserhoilder(A) [m,n]=size(A); Q=eye(n); for i=1:nu=house(A(i:m,i));P=eye(m-i+1)-2/(norm(u)^2)*u*u'; A(i:m,i:n)=P*A(i:m,i:n);PP=blkdiag(eye(i-1),P);Q=Q*PP;endR=triu(A);程序house.mfunction u=house(Ai)Ai(1)=Ai(1)+sign(Ai(1))*norm(Ai);u=Ai/norm(Ai);四、数值实验结果:为了保证在一定阶数下,所做的数值试验具有代表性,我们随机产生3个n阶矩阵,随着n的不断变化,得到如下一系列图:(1)A=rand(4,4)A =0.3027 0.7894 0.1793 0.39220.2364 0.7134 0.4350 0.39360.4690 0.0742 0.9747 0.11150.5381 0.6768 0.7245 0.8650[Q,R]=qr_Houserhoilder(A)Q =-0.3735 -0.5369 0.4602 0.6003-0.2916 -0.5439 -0.7840 -0.0668-0.5786 0.6445 -0.2679 0.4219-0.6639 -0.0208 0.3190 -0.6761R =-0.8106 -0.9951 -1.2387 -0.90010 -0.7781 0.2803 -0.37070 0 -0.2886 0.11800 0 0 -0.3286test_qr3err =1.0e-10 *0.1950 0.0120 0.0119 0.3029-111 1.52 2.53 3.54(2)A=rand(5,5)A =0.7575 0.7005 0.6420 0.8233 0.88770.6373 0.4698 0.5815 0.5290 0.37700.6643 0.3036 0.9458 0.7210 0.12950.1056 0.0254 0.3903 0.1048 0.30270.2813 0.3787 0.5302 0.0627 0.1156[Q,R]=qr_Houserhoilder(A)Q =-0.6161 0.4378 0.2968 -0.5805 -0.0598-0.5183 -0.0381 0.2932 0.7098 -0.3743-0.5403 -0.6729 -0.2440 -0.1030 0.4303-0.0859 -0.1847 -0.5513 -0.2507 -0.7693-0.2288 0.5657 -0.6801 0.2928 0.2817R =-1.2295 -0.9280 -1.3629 -1.1944 -0.86480 0.2940 -0.1496 -0.1288 0.29660 0 -0.4455 0.1231 0.09690 0 0 -0.1847 -0.30310 0 0 0 -0.3388test_qr3err =1.0e-10 *0.1722 0.0119 0.0120 0.4239-111 1.52 2.53 3.54(3)A=rand(6,6)A =0.5961 0.8422 0.4279 0.3473 0.2731 0.34860.1335 0.2860 0.2157 0.6181 0.6306 0.19300.4479 0.0467 0.3284 0.4588 0.4228 0.95480.9430 0.7296 0.7951 0.0920 0.0322 0.31140.7581 0.7811 0.4338 0.1512 0.6185 0.12710.8647 0.1359 0.1735 0.9036 0.4393 0.8964 [Q,R]=qr_Houserhoilder(A)Q =-0.3571 -0.5571 0.2402 0.2038 0.4576 0.5035-0.0800 -0.2530 -0.1880 0.8253 -0.1645 -0.4314-0.2684 0.3532 -0.4622 0.2419 -0.3487 0.6399-0.5650 -0.0874 -0.6012 -0.3372 0.2898 -0.3376-0.4542 -0.3266 0.3151 -0.2560 -0.7185 -0.0773-0.5181 0.6217 0.4823 0.1991 0.2045 -0.1760R =-1.6690 -1.1737 -0.9944 -0.8854 -0.7882 -1.09430 -0.7594 -0.2803 0.3166 -0.0940 0.58290 0 -0.3472 0.1832 0.1390 -0.10880 0 0 0.8020 0.5966 0.50220 0 0 0 -0.4714 -0.02290 0 0 0 0 0.4305test_qr3err =1.0e-10 *0.1518 0.0119 0.0120 0.2845-111 1.52 2.53 3.54五、总结与分析由上面的数值试验结果和绘图可知:比较之下,我自己写的houserhoilder程序误差是最大的,CGS与MGS的误差相差不大,系统自带的matlab软件速度是最快的,我自己写的houserhoilder程序还是最慢的。

数值代数实验报告(1)绝对经典

数值代数实验报告(1)绝对经典

数值代数实验报告Numerical Linear Algebra And ItsApplications学生所在学院:理学院学生所在班级:计算数学10-1学生姓名:戈东潮指导教师:于春肖教务处2012年12月实验一实验名称: Poisson 方程边值问题的五点差分格式 实验时间: 2012年12月13日 星期四 实验成绩: 一、实验目的:通过上机利用Matlab 数学软件实现Poisson 方程的边值问题的五点差分格式的线性代数方程组来认真解读Poisson 方程边值问题的具体思想与方法,使我们掌握得更加深刻,并且做到学习与使用并存,增加学习的实际动手性,不再让学习局限于书本和纸上,而是利用计算机学习来增加我们的学习兴趣。

二、实验内容:利用Poisson 方程来解下列问题:⎪⎩⎪⎨⎧∂∈=∈=∂∂+∂∂-ΩΩ),(,),(),(,sin sin )(y x y x u y x y x y u x u 0222222πππ 其中}{的边界是ΩΩ∂<<∈Ω,1,0),(y x y x 。

边值问题的解释是y x y x u ππsin sin ),(=(1)用例题中模型问题的方法取N=5,10(也可以取N=20)列出五点差分格式的线性代数方程组三、实验过程:系数矩阵A 的Matlab 实现函数如下:%文件名:poisson.m function T=poisson(N) for i=1:N a(i)=4; end b=diag(a); for j=1:N-1b(j,j+1)=-1;endfor j=2:Nb(j,j-1)=-1;endfor i=1:NI(i)=-1;endI=diag(I);for i=1:N:N*Nfor j=1:NT(i+j-1,i:i+N-1)=b(j,:);endendfor i=1:N:N*N-Nfor j=1:NT(i+j-1,i+N:i+N-1+N)=I(j,:);endendfor i=1+N:N:N*Nfor j=1:NT(i+j-1,i-N:i-1)=I(j,:);endend求解常数项b的Matlab实现函数如下:%函数名:constant.mfunction b=constant(N)n=N+1;h=1/n;i=1;for y=1/n:1/n:N/nfor x=1/n:1/n:N/nf(i)=2*pi*pi*sin(pi*x)*sin(pi*y);i=i+1;endendb=h*h*f;b=b';四、实验结果(总结/方案)在Matlab运行窗口输入>>A=poisson(5) Enter输出结果A =Columns 1 through 114 -1 0 0 0 -1 0 0 0 0 0-1 4 -1 0 0 0 -1 0 0 0 00 -1 4 -1 0 0 0 -1 0 0 00 0 -1 4 -1 0 0 0 -1 0 00 0 0 -1 4 0 0 0 0 -1 0-1 0 0 0 0 4 -1 0 0 0 -10 -1 0 0 0 -1 4 -1 0 0 00 0 -1 0 0 0 -1 4 -1 0 00 0 0 -1 0 0 0 -1 4 -1 00 0 0 0 -1 0 0 0 -1 4 00 0 0 0 0 -1 0 0 0 0 40 0 0 0 0 0 -1 0 0 0 -10 0 0 0 0 0 0 -1 0 0 00 0 0 0 0 0 0 0 -1 0 00 0 0 0 0 0 0 0 0 -1 00 0 0 0 0 0 0 0 0 0 -10 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 0 Columns 12 through 220 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 0 0 00 -1 0 0 0 0 0 0 0 0 00 0 -1 0 0 0 0 0 0 0 00 0 0 -1 0 0 0 0 0 0 0 -1 0 0 0 -1 0 0 0 0 0 04 -1 0 0 0 -1 0 0 0 0 0 -1 4 -1 0 0 0 -1 0 0 0 00 -1 4 -1 0 0 0 -1 0 0 00 0 -1 4 0 0 0 0 -1 0 00 0 0 0 4 -1 0 0 0 -1 0 -1 0 0 0 -1 4 -1 0 0 0 -10 -1 0 0 0 -1 4 -1 0 0 00 0 -1 0 0 0 -1 4 -1 0 00 0 0 -1 0 0 0 -1 4 0 00 0 0 0 -1 0 0 0 0 4 -10 0 0 0 0 -1 0 0 0 -1 40 0 0 0 0 0 -1 0 0 0 -10 0 0 0 0 0 0 -1 0 0 00 0 0 0 0 0 0 0 -1 0 0 Columns 23 through 250 0 00 0 00 0 00 0 00 0 00 0 00 0 00 0 00 0 00 0 00 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 -1 0 0 0 -1 0 0 0 -1 0 0 4 -1 0 -1 4 -10 -1 4 在运行窗口输入 >>b=constant(5) Enter 常数项的输出结果b =[ 0.1371 0.2374 0.2742 0.2374 0.1371 0.2374 0.4112 0.4749 0.4112 0.2374 0.2742 0.4749 0.5483 0.4749 0.2742 0.2374 0.4112 0.4749 0.4112 0.2374 0.1371 0.2374 0.2742 0.2374 0.1371]T所以A*u=b 的五点差分格式为j i j i j i j i j i j i f h u u u u u ,,,,,,211114=-----+-+ (i,j=1,2……)实验二实验名称: 用Jacobi 迭代法和SOR 迭代法求解方程组 实验时间: 2012年12月13日 星期四 实验成绩: 一、实验目的:通过上机利用Matlab 数学软件采用Jacobi 迭代法和SOR 迭代法来求解方程组,求解过程中了解两种迭代方法的基本思想与迭代过程,并分析两种迭代方法的收敛性和实用用以及他们的异同点。

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

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

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

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

数据如下:(1) 0.4 0.55 0.65 0.80 0.95 1.05 0.41075 0.57815 0.69675 0.90 1.00 1.25382 求五次Lagrange多项式,和分段三次插值多项式,计算, 的值。

(提示:结果为, )(2) 1 2 3 4 5 6 7 0.368 0.135 0.050 0.018 0.007 0.002 0.001 试构造Lagrange多项式,计算的,值。

(提示:结果为, )二、要求 1、利用Lagrange插值公式编写出插值多项式程序;2、给出插值多项式或分段三次插值多项式的表达式;3、根据节点选取原则,对问题(2)用三点插值或二点插值,其结果如何;4、对此插值问题用Newton插值多项式其结果如何。

Newton 插值多项式如下:其中:三、目的和意义 1、学会常用的插值方法,求函数的近似表达式,以解决其它实际问题;2、明确插值多项式和分段插值多项式各自的优缺点;3、熟悉插值方法的程序编制;4、如果绘出插值函数的曲线,观察其光滑性。

四、实验步骤(1) 0.4 0.55 0.65 0.80 0.951.05 0.41075 0.57815 0.69675 0.90 1.00 1.25382 求五次Lagrange多项式,和分段三次插值多项式,计算, 的值。

数值分析课程的实验报告

数值分析课程的实验报告

实验报告实验课程:学生姓名:学号:专业班级:2012年 6 月 8 日目录---(1)用C语言或C++编程显示字母T与5 (3)---(2)原子弹爆炸能量估计 (11)---(3)城市水管埋多深 (15)---(4)实现PageRank算法 (18)计算机系数值分析实验报告---(1)用C 语言或C ++编程显示字母T 与5学生姓名:学号:专业班级:实验类型:■ 验证 □ 综合 □ 设计□ 创新 实验日期:2012/4/20实验成绩:一、实验目的用C 语言或C ++编程显示字母T 与5。

二、实验基本原理和内容Bezier 曲线生成: 1、确定曲线的阶次;2、计算Bernstein 基函数的表达式:B0,3(t)﹦(1-t)3 ;B1,3(t)﹦3t(1-t)2 ;B2,3(t)﹦3t2(1-t) ;B3,3(t)﹦t3 3、把Bezier 曲线中的Pk 写成分量坐标的形式4、确定一合适的步长;控制t 从0到1变化,求出一系列(x,y)坐标点;将其用小线段顺序连接起来。

三、算法分析对于二维平面的情况,只有x,y 坐标分量,可以给出四点三次Bezier 曲线如下的算法描述: 输入:阶次,3; 控制顶点:4个,(x0,y0),…,(x3,y3) begin x=x0 y=y0 moveto (x,y)for t ﹦0 to 1 step ∆tx ﹦B0,3(t)x0﹢B1,3(t)x1﹢B2,3(t)x2﹢B3,3(t)x3 y ﹦B0,3(t)y0﹢B1,3(t)y1﹢B2,3(t)y2﹢B3,3(t)y3 lineto (x,y) endfor end三次Bezier 曲线例子:设在平面上给定的7个控制点坐标分别为:A (100,300),B (120,200),C (220,200),D (270,100),E (370,100),F (420,200),G (420,300)。

画出其曲线。

n k t t t k n k n t t C t B k n k kn k k nn k ,,1,0]1,0[)1()!(!!)1()(, =∈--=-=--四、实验内容和过程代码:// xx.cpp : Defines the entry point for the application.//#include "stdafx.h"#include "resource.h"#define MAX_LOADSTRING 100// Global Variables:HINSTANCE hInst; // current instanceTCHAR szTitle[MAX_LOADSTRING]; // The title bar text TCHAR szWindowClass[MAX_LOADSTRING]; // The title bar text// Foward declarations of functions included in this code module:ATOM MyRegisterClass(HINSTANCE hInstance);BOOL InitInstance(HINSTANCE, int);LRESULT CALLBACK WndProc(HWND, UINT, WPARAM, LPARAM);LRESULT CALLBACK About(HWND, UINT, WPARAM, LPARAM);//===================================================//===========要我们加的数据================intxT[15][4]={237,237,237,237,237,237,226,143,143,143,143,143,143,143,435,435,435,435,435,435,435,353,339,339,339,339,339,339,339,507,529,552,552,552,576,576, 576,576,570,570,570,570,6,6,6,6,0,0,0,0,24,24,24,48,71,183,183,183,237,237}; intyT[15][4]={620,620,120,120,120,35,24,19,19,19,0,0,0,0,0,0,0,0,19,19,19,23,36,10 9,109,108,620,620,620,620,602,492,492,492,492,492,492,492,662,662,662,662,662,6 62,662,662,492,492,492,492,492,492,492,602,620,620,620,620,620,620};intx5[21][4]={149,149,149,345,345,361,356,368,368,406,368,406,406,397,406,397,397, 382,372,351,351,351,351,142,142,33,142,33,33,32,32,32,32,32,35,44,44,74,109,149 ,149,269,324,324,324,324,264,185,185,165,149,119,119,86,65,42,42,14,0,0,0,0,46, 121,121,205,282,333,333,378,399,399,399,399,381,333,333,288,232,112,112,112,149 ,149};inty5[21][4]={597,597,597,597,597,597,599,606,606,695,606,695,695,702,695,702,702, 681,676,676,676,676,676,676,676,439,676,439,439,438,436,434,434,428,426,426,426 ,426,420,408,408,372,310,208,208,112,37,37,37,37,44,66,66,90,99,99,99,99,87,62, 62,24,0,0,0,0,27,78,78,123,180,256,256,327,372,422,422,468,491,512,512,512,597, 597};//====================================================int APIENTRY WinMain(HINSTANCE hInstance,HINSTANCE hPrevInstance,LPSTR lpCmdLine,intnCmdShow){// TODO: Place code here.MSG msg;HACCEL hAccelTable;// Initialize global stringsLoadString(hInstance, IDS_APP_TITLE, szTitle, MAX_LOADSTRING);LoadString(hInstance, IDC_XX, szWindowClass, MAX_LOADSTRING);MyRegisterClass(hInstance);// Perform application initialization:if (!InitInstance (hInstance, nCmdShow)){return FALSE;}hAccelTable = LoadAccelerators(hInstance, (LPCTSTR)IDC_XX);// Main message loop:while (GetMessage(&msg, NULL, 0, 0)){if (!TranslateAccelerator(msg.hwnd, hAccelTable, &msg)){TranslateMessage(&msg);DispatchMessage(&msg);}}returnmsg.wParam;}//// FUNCTION: MyRegisterClass()//// PURPOSE: Registers the window class.//// COMMENTS://// This function and its usage is only necessary if you want this code// to be compatible with Win32 systems prior to the 'RegisterClassEx'// function that was added to Windows 95. It is important to call this function // so that the application will get 'well formed' small icons associated// with it.//ATOM MyRegisterClass(HINSTANCE hInstance){WNDCLASSEX wcex;wcex.cbSize = sizeof(WNDCLASSEX);wcex.style = CS_HREDRAW | CS_VREDRAW;wcex.lpfnWndProc = (WNDPROC)WndProc;wcex.cbClsExtra = 0;wcex.cbWndExtra = 0;wcex.hInstance = hInstance;wcex.hIcon = LoadIcon(hInstance, (LPCTSTR)IDI_XX);wcex.hCursor = LoadCursor(NULL, IDC_ARROW);wcex.hbrBackground = (HBRUSH)(COLOR_WINDOW+1);wcex.lpszMenuName = (LPCSTR)IDC_XX;wcex.lpszClassName = szWindowClass;wcex.hIconSm = LoadIcon(wcex.hInstance, (LPCTSTR)IDI_SMALL);returnRegisterClassEx(&wcex);}//// FUNCTION: InitInstance(HANDLE, int)//// PURPOSE: Saves instance handle and creates main window//// COMMENTS://// In this function, we save the instance handle in a global variable and // create and display the main program window.//BOOL InitInstance(HINSTANCE hInstance, intnCmdShow){HWND hWnd;hInst = hInstance; // Store instance handle in our global variablehWnd = CreateWindow(szWindowClass, szTitle, WS_OVERLAPPEDWINDOW,CW_USEDEFAULT, 0, CW_USEDEFAULT, 0, NULL, NULL, hInstance, NULL);if (!hWnd){return FALSE;}ShowWindow(hWnd, nCmdShow);UpdateWindow(hWnd);return TRUE;}//// FUNCTION: WndProc(HWND, unsigned, WORD, LONG)//// PURPOSE: Processes messages for the main window.//// WM_COMMAND - process the application menu// WM_PAINT - Paint the main window// WM_DESTROY - post a quit message and return////LRESULT CALLBACK WndProc(HWND hWnd, UINT message, WPARAM wParam, LPARAM lParam){intwmId, wmEvent;PAINTSTRUCT ps;HDC hdc;TCHAR szHello[MAX_LOADSTRING];LoadString(hInst, IDS_HELLO, szHello, MAX_LOADSTRING);switch (message){case WM_COMMAND:wmId = LOWORD(wParam);wmEvent = HIWORD(wParam);// Parse the menu selections:switch (wmId){case IDM_ABOUT:DialogBox(hInst, (LPCTSTR)IDD_ABOUTBOX, hWnd, (DLGPROC)About);break;case IDM_EXIT:DestroyWindow(hWnd);break;default:returnDefWindowProc(hWnd, message, wParam, lParam);}break;case WM_PAINT:hdc = BeginPaint(hWnd, &ps);// TODO: Add any drawing code here...RECT rt;GetClientRect(hWnd, &rt);//=================================================//=======要我们加的代码==========================intx,y,i;double t;for(i=0;i<15;i++){for(t=0.00;t<1;t=t+0.00001){x=(int)(xT[i][0]*(1-t)*(1-t)*(1-t)+xT[i][1]*3*(1-t)*(1-t)*t+xT[i][2]*3*(1-t )*t*t+xT[i][3]*t*t*t)/3;y=300-(int)(yT[i][0]*(1-t)*(1-t)*(1-t)+yT[i][1]*3*(1-t)*(1-t)*t+yT[i][2]*3* (1-t)*t*t+yT[i][3]*t*t*t)/3;SetPixel(hdc, 100+x, y,RGB(255,0,0));}}for(i=0;i<21;i++){for(t=0.00;t<1.0;t=t+0.00001){x=(int)(x5[i][0]*(1-t)*(1-t)*(1-t)+x5[i][1]*3*(1-t)*(1-t)*t+x5[i][2]*3*(1-t )*t*t+x5[i][3]*t*t*t)/4;y=300-(int)(y5[i][0]*(1-t)*(1-t)*(1-t)+y5[i][1]*3*(1-t)*(1-t)*t+y5[i][2]*3* (1-t)*t*t+y5[i][3]*t*t*t)/4;SetPixel(hdc, 300+x, y,RGB(255,0,0));}}//===============================================EndPaint(hWnd, &ps);break;case WM_DESTROY:PostQuitMessage(0);break;default:returnDefWindowProc(hWnd, message, wParam, lParam);}return 0;}// Mesage handler for about box.LRESULT CALLBACK About(HWND hDlg, UINT message, WPARAM wParam, LPARAM lParam) {switch (message){case WM_INITDIALOG:return TRUE;case WM_COMMAND:if (LOWORD(wParam) == IDOK || LOWORD(wParam) == IDCANCEL){EndDialog(hDlg, LOWORD(wParam));return TRUE;}break;}return FALSE;}五、实验结果截图:六、实验心得Bezier曲线的形状是通过一组多边折线(特征多边形)的各顶点唯一地定义出来的,其中只有第一个点和最后一个点在曲线上。

上海大学实验报告模板

上海大学实验报告模板

上海大学实验报告模板一、引言本实验报告旨在介绍上海大学实验报告的模板,帮助学生更好地撰写实验报告。

实验报告是一种科技文献,通过记录实验过程和结果,以及对实验数据的分析和讨论,来总结和分享实验研究的成果。

二、实验报告的结构1.标题:实验报告的标题应简明扼要地概括实验的内容和目的。

2.引言:介绍实验的背景、目的和意义,引起读者的兴趣,并阐明实验的问题和假设。

3.实验方法:详细描述实验的步骤和操作过程,包括实验所使用的设备、材料和方法。

4.实验结果与分析:展示实验的结果数据,并通过图表或数学模型进行分析和解释。

分析的过程应该清晰、逻辑严谨。

5.讨论与结论:对实验结果进行讨论和总结,回答实验的问题和验证假设。

同时提出实验中存在的问题和改进的可能性。

6.参考文献:引用实验中使用的参考文献,包括已发表的科技文献或相关实验报告。

三、实验报告的要求1.格式:实验报告应采用A4纸,页边距为上下左右各2.5厘米,使用12号宋体字体,行间距为1.5倍。

2.内容:实验报告应包含必要的信息,包括实验的背景、目的、方法、结果和分析、讨论和结论等。

3.语言:实验报告应使用准确、简明扼要的语言进行描述,避免使用模糊或不确定的表达方式。

4.图表:实验结果可以通过图表的方式进行展示和分析,但在本实验报告模板中,不允许使用图片。

5.引用:实验报告中引用的参考文献应注明出处,包括作者、标题、出版日期等信息。

引用格式应符合学校的要求。

四、实验报告的撰写步骤1.确定实验的背景和目的:在撰写实验报告之前,需要对实验的背景和目的进行充分的了解和分析,确保能够准确地描述实验的意义和目标。

2.收集实验数据:进行实验过程中,需要仔细记录实验数据和观察结果。

确保实验数据的准确性和完整性。

3.分析实验数据:通过对实验数据的分析和处理,得出实验结论,并将分析过程详细记录在实验报告中。

4.撰写实验报告:根据上述实验报告的结构和要求,逐步撰写实验报告的各个部分,确保内容的完整、准确和逻辑性。

上海大学课程实验报告-数值代数

上海大学课程实验报告-数值代数

课程实验报告COURSE PAPER课程名称:数值代数与计算方法课程号:08305114授课教师:学号:姓名:所属:计算机科学与工程打印时间:2015评语:题目一:算法:1、对于问题一:2、对于问题二:直接编写递归函数程序,算出三个差分方程的10个近似值程序:1)%main.mclc;clear all;a=1;b=-2;c=-3;[x1,x2]=roots(a,b,c)%roots.mfunction [x1,x2]=roots(a,b,c)d=sqrt(b*b-4*a*c);if d>0x1=(-2*c)/(b+d);x2=(-b-d)/(2*a);elseif d<0x1=(-b+d)/(2*a);x2=(-2*c)/(b-d);end2)%solu.mfunction [X,R,P,Q]=solu(X0,R0,P0,P1,Q0,Q1) X(1)=X0;R(1)=R0;P(1)=P0;P(2)=P1;Q(1)=Q0;Q(2)=Q1;for i=1:9X(i+1)=X(i)/2;X(i)=X(i+1);endfor i=1:9R(i+1)=R(i)/2;R(i)=R(i+1);endfor i=3:10P(i)=3/2*P(i-1)-1/2*P(i-2);P(i-1)=P(i);P(i-2)=P(i-1);endfor i=3:10Q(i)=5/2*Q(i-1)-Q(i-2);Q(i-1)=Q(i);Q(i-2)=Q(i-1);end%Xn.mclc;clear all;X0=1;R0=0.994;P0=1;P1=0.497;Q0=1;Q1=0.497;[X,R,P,Q]=solu(X0,R0,P0,P1,Q0,Q1);for i=1:10x(i)=i;endplot(x,X-R,'r*');hold on结果:1)x1 =3x2 =-12)题目二:1、此题我们以第一个式子为例,求出g(Xn+1)= (2/(2+3*Xn-Xn.^3))^(1/2)2、利用二分法:3、利用牛顿迭代法,并根据题目提示的立方根算法:此题我们以第三个式子为例:,设A=71/3 ,求出x的3次方便是近似值。

上海大学实验报告

上海大学实验报告

竭诚为您提供优质文档/双击可除上海大学实验报告篇一:上海大学操作系统(二)实验报告(全)评分:shAnghAIunIVeRsITY操作系统实验报告学院计算机工程与科学专业计算机科学与技术学号学生姓名《计算机操作系统》实验一报告实验一题目:操作系统的进程调度姓名:张佳慧学号:12122544实验日期:20XX.1实验环境:microsoftVisualstudio实验目的:进程是操作系统最重要的概念之一,进程调度又是操作系统核心的主要内容。

本实习要求学生独立地用高级语言编写和调试一个简单的进程调度程序。

调度算法可任意选择或自行设计。

例如,简单轮转法和优先数法等。

本实习可加深对于进程调度和各种调度算法的理解。

实验内容:1、设计一个有n个进程工行的进程调度程序。

每个进程由一个进程控制块(pcb)表示。

进程控制块通常应包含下述信息:进程名、进程优先数、进程需要运行的时间、占用cpu的时间以及进程的状态等,且可按调度算法的不同而增删。

2、调度程序应包含2~3种不同的调度算法,运行时可任意选一种,以利于各种算法的分析比较。

3、系统应能显示或打印各进程状态和参数的变化情况,便于观察诸进程的调度过程。

操作过程:1、本程序可选用优先数法或简单轮转法对五个进程进行调度。

每个进程处于运行R(run)、就绪w(wait)和完成F(finish)三种状态之一,并假设起始状态都是就绪状态w。

为了便于处理,程序进程的运行时间以时间片为单位计算。

进程控制块结构如下:进程控制块结构如下:pcb进程标识数链指针优先数/轮转时间片数占用cpu时间片数进程所需时间片数进程状态进程控制块链结构如下:其中:Run—当前运行进程指针;heAD—进程就绪链链首指针;TAID—进程就绪链链尾指针。

2、算法与框图(1)优先数法。

进程就绪链按优先数大小从高到低排列,链首进程首先投入运行。

每过一个时间片,运行进程所需运行的时间片数减1,说明它已运行了一个时间片,优先数也减3,理由是该进程如果在一个时间片中完成不了,优先级应该降低一级。

高等数学数学实验报告(两篇)2024

高等数学数学实验报告(两篇)2024

引言概述:高等数学数学实验报告(二)旨在对高等数学的相关实验进行探究与研究。

本次实验报告共分为五个大点,每个大点讨论了不同的实验内容。

在每个大点下,我们进一步细分了五到九个小点,对实验过程、数据收集、数据分析等进行了详细描述。

通过本次实验,我们可以更好地理解高等数学的概念和应用。

正文内容:一、微分方程实验1.利用欧拉法求解微分方程a.介绍欧拉法的原理和步骤b.详细阐述欧拉法在实际问题中的应用c.给出具体的实例,展示欧拉法的计算步骤2.应用微分方程建立模型求解实际问题a.介绍微分方程模型的建立方法b.给出一个具体的实际问题,使用微分方程建立模型c.详细阐述模型求解步骤和结果分析3.使用MATLAB求解微分方程a.MATLAB求解微分方程的基本语法和函数b.给出一个具体的微分方程问题,在MATLAB中进行求解c.分析结果的准确性和稳定性二、级数实验1.了解级数的概念和性质a.简要介绍级数的定义和基本概念b.阐述级数收敛和发散的判别法c.讨论级数的性质和重要定理2.使用级数展开函数a.介绍级数展开函数的原理和步骤b.给出一个函数,使用级数展开进行近似计算c.分析级数近似计算的精确度和效果3.级数的收敛性与运算a.讨论级数收敛性的判别法b.介绍级数的运算性质和求和法则c.给出具体的例题,进行级数的运算和求和三、多元函数极值与最值实验1.多元函数的极值点求解a.介绍多元函数的极值点的定义和求解方法b.给出一个多元函数的实例,详细阐述求解过程c.分析极值点对应的函数值和意义2.多元函数的条件极值与最值a.讨论多元函数的条件极值的判定法b.给出一个具体的多元函数,求解其条件极值和最值c.分析条件极值和最值对应的函数值和意义3.利用MATLAB进行多元函数极值与最值的计算a.MATLAB求解多元函数极值与最值的基本语法和函数b.给出一个多元函数的具体问题,在MATLAB中进行求解c.分析结果的准确性和可行性四、曲线积分与曲面积分实验1.曲线积分的计算方法与应用a.介绍曲线积分的定义和计算方法b.给出一个具体的曲线积分问题,详细阐述计算过程c.分析曲线积分结果的几何意义2.曲线积分的应用举例a.讨论曲线积分在实际问题中的应用b.给出一个实际问题,使用曲线积分进行求解c.分析曲线积分结果的实际意义和应用价值3.曲面积分的计算方法与应用a.介绍曲面积分的定义和计算方法b.给出一个具体的曲面积分问题,详细阐述计算过程c.分析曲面积分结果的几何意义五、空间解析几何实验1.空间曲线的参数方程表示与性质a.介绍空间曲线的参数方程表示和性质b.给出一个具体的空间曲线,转化为参数方程表示c.分析参数方程对应的几何意义和性质2.平面与空间直线的位置关系a.讨论平面与空间直线的位置关系的判定方法b.给出一个具体的平面与空间直线的问题,判定其位置关系c.分析位置关系对应的几何意义和应用实例3.空间直线与平面的夹角和距离计算a.介绍空间直线与平面的夹角和距离的计算方法b.给出一个具体的空间直线和平面,计算其夹角和距离c.分析夹角和距离计算结果的几何意义总结:通过本次高等数学数学实验报告(二),我们深入了解了微分方程、级数、多元函数极值与最值、曲线积分、曲面积分以及空间解析几何的相关概念和应用。

数值代数实验报告1

数值代数实验报告1

暨南大学本科实验报告专用纸课程名称数值代数成绩评定实验项目名称(填写实验所属章节名称) 指导教师刘娟实验项目编号0701******* 实验项目类型验证实验地点机房学生姓名王伟文学号2010051727学院信息科学技术学院数学系信息与计算科学专业2010级实验时间2011年9月 1 日~12月30日温度24℃第一章1.实验选题:(写出你在该项目所选的实验题目)第一章上机练习12.谈谈你对该算法的理解:(简单谈一下你是如何理解该算法的?)对算法的理解:先将84阶的矩阵A分解为一个下三角矩阵L和上三角矩阵U,先考虑下三角形方程组Ly=b,利用前代法解出y,再考虑上三角形方程组Ux=y;利用回代法解出x。

列主元gauss消去法Ax b PA LU Ly Pb Ux y=⇔===,,;3.实验内容(将实验程序及其实验结果粘贴,最好对程序各部分注释清楚,比如设置了哪些函数,这些函数的输入输出是什么,具有什么功能?)实验程序对矩阵A进行LU分解的程序function [ L,U ] = LUfac( A )for k=1:n-1A(k+1:n,k)=A(k+1:n,k)/A(k,k);A(k+1:n,k+1:n)=A(k+1:n,k+1:n)-A(k+1:n,k)*A(k,k+1:n);endL=tril(A,0);for i=1:nL(i,i)=1;endU=triu(A,0);end利用前代法解出y值的程序function [ b ] = TSL( L,b )n=size(L,1);for j=1:n-1b(j)=b(j)/L(j,j);b(j+1:n)=b(j+1:n)-b(j)*L(j+1:n,j);endb(n)=b(n)/L(n,n);end利用回代法解出x值的程序function [ b ] = TSU( U,b )n=size(U,1);for j=n:-1:2b(j)=b(j)/U(j,j);b(1:j-1)=b(1:j-1)-b(j)*U(1:j-1,j);endb(1)=b(1)/U(1,1);end主函数程序(‘生成84阶的矩阵A’)A=eye(84);A=6*A;for i=2:84A(i,i-1)=8;A(i-1,i)=1;end(‘生成84乘1的矩阵b’)b=ones(84,1);b=b*15;b(1)=7;b(84)=14;[L,U]=LUfac(A);(‘调用函数LUfac对矩阵A进行分解’)y=TSL(L,b);(‘调用函数TSL求解Ly=b方程’)x=TSU(U,y);(‘调用函数TSU求解UX=b方程’)用MATLAB求解得到的结果x’ans =1.0e+008 *Columns 1 through 70.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 8 through 140.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 15 through 210.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 22 through 280.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 29 through 350.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 36 through 420.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 43 through 490.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 50 through 560.0000 0.0000 0.0000 0.0000 0.0000 -0.0000 0.0000 Columns 57 through 63-0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 Columns 64 through 700.0000 -0.0000 0.0000 -0.0001 0.0002 -0.0003 0.0007 Columns 71 through 77-0.0013 0.0026 -0.0052 0.0105 -0.0209 0.0419 -0.0836 Columns 78 through 840.1665 -0.3303 0.6501 -1.2582 2.3487 -4.0263 5.3684列主元gauss消去法function [L,U,P]=Lufac(A)n=size(A,1);P=eye(n,n);for i=1:n-1[r,m]=max(abs(A(i:n,i)));m=m+i-1;A([i,m],:)=A([m,i],:);P([i,m])=P([m,i]);if A(i,i)~=0A(i+1:n,i)=A(i+1:n,i)/ A(i,i);A(i+1:n,i+1:n)= A(i+1:n,i+1:n)-A(i+1:n,i)*A(i,i+1:n); endendU=triu(A);L=tril(A,-1)+eye(n,n);endA=eye(84);A=6*A;for i=2:84A(i,i-1)=8; A(i-1,i)=1; endb=ones(84,1);b=b*15;b(1)=7;b(84)=14;[L,U,P]=Lufac(A); b=P*b;y=TSL(L,b);x=TSU(U,y);求解结果为x =1.0e+025 *0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00010.0002-0.00040.0007-0.00150.0030-0.00600.0119-0.02380.0475-0.09460.1878-0.36960.7154-1.33552.2894-3.05251.教师评语、评分:(请认真对待每次数值代数实验,期末交实验报告,将计算实验成绩)2.实验选题:第三章上机习题1(3)。

大学数学实验报告总结(3篇)

大学数学实验报告总结(3篇)

第1篇一、实验背景随着科学技术的不断发展,数学在各个领域的应用日益广泛。

为了提高学生运用数学知识解决实际问题的能力,本实验课程旨在通过一系列数学实验,让学生深入理解数学理论,掌握数学软件的使用,并培养创新思维和团队协作精神。

二、实验目的1. 深入理解数学理论知识,提高数学应用能力。

2. 掌握数学软件(如MATLAB、Mathematica等)的基本操作和编程技巧。

3. 培养创新思维和团队协作精神,提高实践能力。

4. 通过实验,验证数学理论在实际问题中的应用价值。

三、实验内容本实验课程共分为以下几个部分:1. 数值分析实验:包括数值微分、数值积分、线性方程组的求解等。

2. 线性代数实验:包括矩阵运算、特征值与特征向量、线性方程组的求解等。

3. 概率论与数理统计实验:包括随机变量及其分布、参数估计、假设检验等。

4. 运筹学实验:包括线性规划、整数规划、网络流等。

5. 高等数学实验:包括常微分方程、偏微分方程、复变函数等。

四、实验过程1. 实验准备:查阅相关资料,了解实验原理和方法,明确实验目的和步骤。

2. 实验实施:按照实验指导书的要求,利用数学软件进行实验操作,记录实验数据。

3. 数据分析:对实验数据进行处理和分析,验证数学理论在实际问题中的应用。

4. 实验报告撰写:总结实验过程、结果和心得体会,撰写实验报告。

五、实验结果与分析1. 数值分析实验:通过数值微分、数值积分等方法,验证了数值方法在求解实际问题中的有效性。

例如,在求解非线性方程组时,采用了牛顿迭代法,成功找到了方程的近似解。

2. 线性代数实验:通过矩阵运算、特征值与特征向量等方法,解决了实际工程问题中的线性方程组求解问题。

例如,在求解电路分析问题时,利用矩阵方法求得了电路的电压和电流分布。

3. 概率论与数理统计实验:通过随机变量及其分布、参数估计、假设检验等方法,分析了实际问题中的数据,得出了可靠的结论。

例如,在产品质量检测中,利用假设检验方法判断了产品是否合格。

数值实验报告 - 侯瑞祥

数值实验报告 - 侯瑞祥

本科实验报告
课程名称:计算机数值方法
实验项目:方程求根
实验地点:致远楼B301
专业班级:软件1216 学号:2012005403 学生姓名:侯瑞祥
指导教师:田华
2014年 5 月21 日
本科实验报告
课程名称:计算机数值方法
实验项目:线性方程组的直接求法
实验地点:致远楼B301
专业班级:软件1216 学号:2012005403 学生姓名:侯瑞祥
指导教师:田华
2014年 5 月21 日
本科实验报告
课程名称:计算机数值方法
实验项目:线性方程组的迭代求法
实验地点:致远楼B301
专业班级:软件1216 学号:2012005403 学生姓名:侯瑞祥
指导教师:田华
2014年 5 月21 日
分析:使用高斯-赛德尔和雅克比迭代都可以求出方程组的解,但是利用高斯-赛德尔迭代法所需的迭代次数比雅克比迭代少,收敛快,能够更早的达到精度要求,但是
本科实验报告
课程名称:计算机数值方法
实验项目:代数插值和最小二乘法拟合实验地点:致远楼B301
专业班级:软件1216 学号:2012005403 学生姓名:侯瑞祥
指导教师:田华
2014年 5 月21 日。

数值代数上机实验报告

数值代数上机实验报告

数值代数课程设计实验报告姓名: 班级: 学号: 实验日期:一、实验名称 代数的数值解法 二、实验环境实验一、平方根法与改良平方根法一、实验要求:用熟悉的运算机语言将不选主元和列主元Gasuss 消元法编写成通用的子程序,然后用编写的程序求解以下方程组⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡--⨯14151515157681681681681681612321n n n n n x x x x x x 用所编的程序别离求解40、84、120阶方程组的解。

二、算法描述及实验步骤GAuss 程序如下:(1)求A 的三角分解:LU A =;(2)求解b y =L 得y ; (3)求解y x =U 得x ;列主元Gasuss 消元法程序如下: 1求A 的列主元分解:LU PA =;2求解b y P L =得y ;3求解y x U 得x ;三、调试进程及实验结果:%----------------方程系数---------------->> A1=Sanduijiaozhen(8,6,1,40); >> A2=Sanduijiaozhen(8,6,1,84); >> A3=Sanduijiaozhen(8,6,1,120); >> b1(1)=7;b2(1)=7;b3(1)=7; >> for i=2:39 b1(i)=15; end>> b1(40)=14; >> for i=2:83 b2(i)=15; end>> b2(40)=14; >> for i=2:119 b1(i)=15; end>> b3(120)=14;%----------------方程解---------------->> x11=GAuss(A1,b1') >> x12=GAuss Zhu (A1,b1') >> x21=GAuss(A2,b2') >> x22=GAuss Zhu (A3,b3') >> x31=GAuss(A3,b3') >> x32=GAuss Zhu (A3,b3') 运行结果:(n=40) GAuss 消元法的解即为x11 =列主元GAuss消元法的解即为x12 =1 1 1 1 1 1 1 1 1 1 111111111111111111111111111111六、源程序:function A=Sanduijiaozhen(a,b,c,n)%生成n阶以a,b,c为元素的三对角阵A=diag(b*ones(1,n),0)+diag(c*ones(1,n-1),1)+diag(a*ones(1,n-1),-1);function x=GAuss(A,b)n=length(b);x=b;%-------分解---------------for i=1:n-1for j=i+1:nmi=A(j,i)/A(i,i);b(j)=b(j)-mi*b(i);for k=i:nA(j,k)=A(j,k)-mi*A(i,k);endAB=[A,b];endend%-----------回代------------------x(n)=b(n)/A(n,n);for i=n-1:-1:1s=0;for j=i+1:ns=s+A(i,j)*x(j);endx(i)=(b(i)-s)/A(i,i);endfunction x=GAussZhu(A,b)n=length(b);x=b;%----------------------选主元---------------------for k=1:n-1a_max=0;for i=k:nif abs(A(i,k))>a_maxa_max=abs(A(i,k));r=i;endendif r>kfor j=k:nz=A(k,j);A(k,j)=A(r,j);A(r,j)=z;endz=b(k);b(k)=b(r);b(r)=z;end%--------------消元-----------------for i=k+1:nm=A(i,k)/A(k,k);for j=k:nA(i,j)=A(i,j)-m*A(k,j);endb(i)=b(i)-m*b(k);endendif abs(A(n,n))==0return;endAbZhu=[A,b];%----------------回代-----------------------x(n)=b(n)/A(n,n);for i=n-1:-1:1for j=i+1:nb(i)=b(i)-A(i,j)*x(j);endx(i)=b(i)/A(i,i);end实验二、平方根法与改良平方根法一、实验要求:用运算机语言将平方根法和改良的平方根法编成通用的子程序,然后用编写的程序求解对称正定方程组100阶方程组AX=b,二、算法描述及实验步骤:平方根法函数程序如下:一、求A 的Cholesky 分解:L L A T=; 二、求解b y =L 得y ;3、求解y x =TL 得x ;改良平方根法函数程序如下:一、求A 的Cholesky 分解:T=LDL A ; 二、求解b y =L 得y ; 3、求解y x =TDL 得x ;三、调试进程及实验结果:clear;clc;%----------------方程系数---------------->> A=Sanduijiaozhen(1,10,1,100); >> b(1)=11; >> for i=2:99 b(i)=12; end>> b(100)=11;>> x1=Cholesky(A,b') >> x2=GJCholesky(A,b')运行结果:平方根法的解即为 x1 =改良平方根法解得的解即为x2 =四、源程序:function x=Cholesky(A,b)n=size(A);n=n(1);% x=A^-1*b;% disp('Matlab自带解即为x');%-----------------Cholesky分解-------------------for k=1:nA(k,k)=sqrt(A(k,k));A(k+1:n,k)=A(k+1:n,k)/A(k,k);for j=k+1:n;A(j:n,j)=A(j:n,j)-A(j:n,k)*A(j,k);endend%------------------前代法求解Ly=b----------------------------for j=1:n-1b(j)=b(j)/A(j,j);b(j+1:n)=b(j+1:n)-b(j)*A(j+1:n,j);endb(n)=b(n)/A(n,n);%-----------------回代法求解L'x=y-----------------------------A=A';for j=n:-1:2b(j)=b(j)/A(j,j);b(1:j-1)=b(1:j-1)-b(j)*A(1:j-1,j);endb(1)=b(1)/A(1,1);disp('平方根法的解即为');function b=GJCholesky(A,b)n=size(A);n=n(1);v=zeros(n,1);%----------------------LDL'分解-----------------------------for j=1:nfor i=1:j-1v(i)=A(j,i)*A(i,i);endA(j,j)=A(j,j)-A(j,1:j-1)*v(1:j-1);A(j+1:n,j)=(A(j+1:n,j)-A(j+1:n,1:j-1)*v(1:j-1))/A(j,j);endB=diag(A);D=zeros(n);for i=1:nD(i,i)=B(i);A(i,i)=1;end%-------------------前代法---------------------------A=tril(A); %取得L和Dfor j=1:n-1b(j)=b(j)/A(j,j);b(j+1:n)=b(j+1:n)-b(j)*A(j+1:n,j);endb(n)=b(n)/A(n,n);%-----------------回代法-----------------------------A=D*(A');for j=n:-1:2b(j)=b(j)/A(j,j);b(1:j-1)=b(1:j-1)-b(j)*A(1:j-1,j);endb(1)=b(1)/A(1,1);disp('改良平方根法解得的解即为');实验三、二次多项式拟合一、实验要求:用运算机语言编制利用QR分解求解线性最小二乘问题的通用子程序,用编写的程序求解一个二次多项式使在残向量的范数最小的意义下拟合下面的数据t-1 0iyi二、算法描述及实验步骤:QR分解求解程序如下:一、求A 的QR 分解; 二、计算b c 11T=Q ;3、求解上三角方程1c x =R 得x ;三、调试进程及实验结果:>> t=[-1 0 ]; >> y=[ ]; >> plot(t,y,'r*');>> legend('实验数据(ti,yi)'); >> xlabel('t'), ylabel('y');>> title('二次多项式拟合的数据点(ti,yi)的散点图');运行后屏幕显示数据的散点图(略).(3)编写以下MATLAB 程序计算)(x f 在),(i i y x 处的函数值,即输入程序 >> syms a b c >> t=[-1 0 ]; >> fi=a.*t.^2+ b.*t+c%运行后屏幕显示关于 ,,a b c 的线性方程组fi =[a-b+c,9/16*a-3/4*b+c,1/4*a-1/2*b+c,c,1/16*a+1/4*b+c,1/4*a+1/2*b+c,9/16*a+3/4*b +c]编写构造残向量2范数的MATLAB 程序>> y=[ ]; >> y=[ ];>> fy=fi-y; fy2=fy.^2; J=sum(fy.^2); 运行后屏幕显示误差平方和如下 J=(a-b+c-1)^2+(9/16*a-3/4*b+c-13/16)^2+(1/4*a-1/2*b+c-3/4)^2+(c-1)^2+(1/16*a+1/4*b+c-21/16)^2+(1/4*a+1/2*b+c-7/4)^2+(9/16*a+3/4*b+c-37/16)^2为求,,a b c 使J 达到最小,只需利用极值的必要条件0J a ∂=∂,0J b ∂=∂,0Jc∂=∂,取得关于,,a b c 的线性方程组,这能够由下面的MATLAB 程序完成,即输入程序 >> Ja1=diff(J,a); Ja2=diff(J,b); Ja3=diff(J,c);>> Ja11=simple(Ja1), Ja21=simple(Ja2), Ja31=simple(Ja3) 运行后屏幕显示J 别离对,,a b c 的偏导数如下 Ja11 =451/128*a-63/32*b+43/8*c-887/128 Ja21 =-63/32*a+43/8*b-3/2*c-61/32Ja31 =43/8*a-3/2*b+14*c-143/8解线性方程组112131000Ja Ja Ja ===,,,输入以下程序 >> A=[451/128, -63/32, -3/2 ;-63/32,43/8,-3/2;43/8,-3/2,14]; >> B=[887/128,61/32,143/8]; >> C=B/A, f=poly2sym(C)运行后屏幕显示拟合函数f 及其系数C 如下 C =f =924/2999*x^2+10301/11996*x+4204/2999 故所求的拟合曲线为2()0.30810.8581 1.4018f x x x =++四、源程序:>> t=[-1 0 ]; >> y=[ ]; >> plot(t,y,'r*');>> legend('实验数据(ti,yi)'); >> xlabel('t'), ylabel('y');>> title('二次多项式拟合的数据点(ti,yi)的散点图'); >> syms a b c>> t=[-1 0 ]; >> fi=a.*t.^2+ b.*t+c fi =[ a-b+c, 9/16*a-3/4*b+c, 1/4*a-1/2*b+c, c, 1/16*a+1/4*b+c, 1/4*a+1/2*b+c, 9/16*a+3/4*b+c]>> y=[ ]; >> y=[ ];>> fy=fi-y; fy2=fy.^2; J=sum(fy.^2) J =(a-b+c-1)^2+(9/16*a-3/4*b+c-13/16)^2+(1/4*a-1/2*b+c-3/4)^2+(c-1)^2+(1/16*a+1/4*b+c-21/16)^2+(1/4*a+1/2*b+c-7/4)^2+(9/16*a+3/4*b+c-37/16)^2>> Ja1=diff(J,a); Ja2=diff(J,b); Ja3=diff(J,c);>> Ja11=simple(Ja1), Ja21=simple(Ja2), Ja31=simple(Ja3)Ja11 =451/128*a-63/32*b+43/8*c-887/128Ja21 =-63/32*a+43/8*b-3/2*c-61/32Ja31 =43/8*a-3/2*b+14*c-143/8>> A=[451/128, -63/32, -3/2 ;-63/32,43/8,-3/2;43/8,-3/2,14]; >> B=[887/128,61/32,143/8];>> C=B/A, f=poly2sym(C)C =f =924/2999*x^2+10301/11996*x+4204/2999>>。

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

课程实验报告COURSE PAPER课程名称:数值代数与计算方法课程号:08305114授课教师:学号:姓名:所属:计算机科学与工程打印时间:2015评语:题目一:算法:1、对于问题一:2、对于问题二:直接编写递归函数程序,算出三个差分方程的10个近似值程序:1)%main.mclc;clear all;a=1;b=-2;c=-3;[x1,x2]=roots(a,b,c)%roots.mfunction [x1,x2]=roots(a,b,c)d=sqrt(b*b-4*a*c);if d>0x1=(-2*c)/(b+d);x2=(-b-d)/(2*a);elseif d<0x1=(-b+d)/(2*a);x2=(-2*c)/(b-d);end2)%solu.mfunction [X,R,P,Q]=solu(X0,R0,P0,P1,Q0,Q1) X(1)=X0;R(1)=R0;P(1)=P0;P(2)=P1;Q(1)=Q0;Q(2)=Q1;for i=1:9X(i+1)=X(i)/2;X(i)=X(i+1);endfor i=1:9R(i+1)=R(i)/2;R(i)=R(i+1);endfor i=3:10P(i)=3/2*P(i-1)-1/2*P(i-2);P(i-1)=P(i);P(i-2)=P(i-1);endfor i=3:10Q(i)=5/2*Q(i-1)-Q(i-2);Q(i-1)=Q(i);Q(i-2)=Q(i-1);end%Xn.mclc;clear all;X0=1;R0=0.994;P0=1;P1=0.497;Q0=1;Q1=0.497;[X,R,P,Q]=solu(X0,R0,P0,P1,Q0,Q1);for i=1:10x(i)=i;endplot(x,X-R,'r*');hold on结果:1)x1 =3x2 =-12)题目二:1、此题我们以第一个式子为例,求出g(Xn+1)= (2/(2+3*Xn-Xn.^3))^(1/2)2、利用二分法:3、利用牛顿迭代法,并根据题目提示的立方根算法:此题我们以第三个式子为例:,设A=71/3 ,求出x的3次方便是近似值。

1、%fixpt.mfunction [k,p,err,P]=fixpt(p0,tol,max1)P(1)=p0;for k=2:max1P(k)=g(P(k-1));err=abs(P(k)-P(k-1));relerr=err/(abs(P(k))+eps);p=P(k);if (err<tol)||(relerr<tol),break;endendif k==max1disp('maximum number of iterations exceeded')endP=P';%g.mfunction y=g(x)y=sqrt(2/(2+3*x-x.^3));%iteration.mclc;clear all;p0=0.5;tol=1.0e-9;max1=20;[k,p,err,P]=fixpt(p0,tol,max1);2、二分法:%bisect.mfunction [c,err,yc,k,an,cn,bn,ycn]=bisect(a,b,delta) ya=f(a);yb=f(b);an(1)=a;bn(1)=b;if ya*yb>0disp('No Root in this district!');endmax1=1+round((log(b-a)-log(delta))/log(2));for k=1:max1c=(a+b)/2;cn(k)=c;yc=f(c);ycn(k)=yc;if yc==0a=c;an(k+1)=a;b=c;bn(k+1)=b;elseif yb*yc>0b=c;bn(k+1)=b;yb=yc;elsea=c;an(k+1)=a;ya=yc;endif b-a<deltabreak;endendc=(a+b)/2;cn(max1)=c;err=abs(b-a);yc=f(c);ycn(max1)=yc;for i=1:9k(i)=i-1;end%f.mfunction h=f(x)h=x*sin(x)-1;%main.mclc;clear all;a=0;b=2;delta=0.000001;[c,err,yc,k,an,cn,bn,ycn]=bisect(a,b,delta);4、牛顿迭代法求解:%newton.mfunction [p0,err,k,y]=newton(A,p0,delta,epsilon,max1)for k=1:max1p1=(2.*p0+A/p0.^2)/3;err=abs(p1-p0);relerr=2*err/(abs(p1)+delta);p0=p1;y=f(p0,A);if(err<delta)||(relerr<delta)||(abs(y)<epsilon),break,end end%f.mfunction y=f(x,A)y=x.^3-A;%main.mclc;clear all;A=(7).^(1/3);p0=2;delta=0.001;epsilon=0.000001;max1=30;[p0,err,k,y]=newton(A,p0,delta,epsilon,max1);x=p0^3;vpa(x,10)结果:1、2、二分法3、即71/3的近似值为1.912931183题目三:算法:▲回代算法:程序:%backsub.mfunction X=backsub(A,B)n=length(B);X=zeros(n,1);X(n)=B(n)/A(n,n);for k=n-1:-1:1X(k)=(B(k)-A(k,k+1:n)*X(k+1:n))/A(k,k); end%main.mclear all;clc;A=[3,-2,1,-1;0,4,-1,2;0,0,2,3;0,0,0,5];B=[8,-3,11,15]';X=backsub(A,B)det(A)结果:X =2-213ans =120即x1=2, x2=-2, x3=1, x4=3;|A|=120.题目四:求解方程组:算法:先将方程组系数矩阵进行上三角变换,再进行回代,最终求出方程组解。

上三角变换的算法:高斯消元法,AX=B,设增广矩阵为[A|B],进行初等行变换,通过选取主元行并消元,构成使增广矩阵的A矩阵构成上三角矩阵。

选主元策略常用偏序选主元策略:首先检查位于主对角线或主对角线下方第p列的所有元素,确定行k,它的元素的绝对值最大,即|a kp|=max{|a pp|,|a p+1 p|,…,|a N-1 p|,|a N p|}如果k>p,则交换第k行和第p行。

程序:%uptrbk.mfunction X=uptrbk(A,B)[N N1]=size(A);X=zeros(N,1);C=zeros(1,N+1);Aug=[A B];for p=1:N-1[Y,j]=max(abs(Aug(p:N,p)));C=Aug(p,:);Aug(p,:)=Aug(j+p-1,:);Aug(j+p-1,:)=C;if Aug(p,p)==0'A was singular. NO unique solution'breakendfor k=p+1:Nm=Aug(k,p)/Aug(p,p);Aug(k,p:N+1)=Aug(k,p:N+1)-m*Aug(p,p:N+1);endendX=backsub(Aug(1:N,1:N),Aug(1:N,N+1));%main.mclear all;clc;A=[2,4,-6;1,5,3;1,3,2];B=[-4,10,5]';X=uptrbk(A,B)det(A)结果:X =-321ans =18即x1=-3, x2= 2, x3=1; |A|=18. 题目五:算法:如下线性方程组:1、雅可比迭代法2、高斯-赛德尔迭代法:在此题中我们分别用三组初始值,P1=[x1=0,x2=0,x3=0];P2=[x1=1,x2=2,x3=3];P3=[x1=5,x2=5,x3=5];分别用雅可比和高斯进行迭代,得出结果并比较。

程序:%jacobi.m//雅可比迭代函数function X=jacobi(A,B,P,delta,max1)N=length(B);for k=1:max1for j=1:NX(j)=(B(j)-A(j,[1:j-1,j+1:N])*P([1:j-1,j+1:N]))/A(j,j); enderr=abs(norm(X'-P));relerr=err/(norm(X)+eps);P=X';if(err<delta)||(relerr<delta)breakendendX=X';%gseid.m//高斯-赛德尔迭代函数function X=gseid(A,B,P,delta,max1)N=length(B);for k=1:max1for j=1:Nif j==1X(1)=(B(1)-A(1,2:N)*P(2:N))/A(1,1);elseif j==NX(N)=(B(N)-A(N,1:N-1)*(X(1:N-1))')/A(N,N);elseX(j)=(B(j)-A(j,1:j-1)*X(1:j-1)'-A(j,j+1:N)*P(j+1:N))/A(j,j); endenderr=abs(norm(X'-P));relerr=err/(norm(X)+eps);P=X';if(err<delta)||(relerr<delta)breakendendX=X';%main.mclear all;clc;delta=0.0000001;max1=100;A=[1,0,1;-1,1,0;1,2,-3];B=[2,0,0]';P1=[0,0,0]';P2=[1,2,3]';P3=[5,5,5]';X1=jacobi(A,B,P1,delta,max1);X1=vpa(X1',10)X2=jacobi(A,B,P2,delta,max1);X2=vpa(X2',10)X3=jacobi(A,B,P3,delta,max1);X3=vpa(X3',10)X4=gseid(A,B,P1,delta,max1);X4=vpa(X4',10)X5=gseid(A,B,P2,delta,max1);X5=vpa(X5',10)X6=gseid(A,B,P3,delta,max1);X6=vpa(X6',10)结果:X1 =[ 1.001606818, 1.004722712, 1.003011523]X2 =[ 1.0022603, 0.9955518438, 0.9943430273]X3 =[ 0.9935727281, 0.9811091525, 0.9879539061]X4 =[ 0, 0, 0]X5 =[ 3.0, 3.0, 3.0]X6 =[ 5.0, 5.0, 5.0]正确结果为:x1=x2=x3=1,发现雅可比迭代对于此方程组有效。

相关文档
最新文档