数值代数上机实验报告
数值分析上机实验报告3
实验报告三题目：函数逼近——曲线拟合目的：掌握曲线拟合基本使用方法数学原理：[P,S]=polyfit(x,y,3)其中x，y为取样值，3为得出的结果的最高次数。
P为对应次数的系数，S为误差值向量，其中x,y是等长的向量，P是一个长度为m+1的向量。
结果分析和讨论：23.观察物体的运动，得出时间t与距离s的关系如表，求运动方程。
t=[0,0.9,1.9,3.0,3.9,5.0];s=[0,10,30,50,80,110];[P,S]=polyfit(t,s,5)P =-0.5432 6.4647 -26.5609 46.1436 -13.2601 -0.0000S =R: [6x6 double]df: 0normr: 1.2579e-012所以得到方程为：-13.2601x46.1436x-26.5609x6.4647x-0.5432x2345++24.在某化学反应堆里，根据实验所得分解物的质量分数y与时间t的关系，用最小拟合求y=F（t）；>> x=0:5:55;y=[0,1.27,2.16,2.86,3.44,3.87,4.15,4.37,4.51,4.58,4.62,4.64];>> [P,S]=polyfit(x,y,5)P =0.0000 -0.0000 0.0002 -0.0084 0.2851 0.0082S =R: [6x6 double]df: 6normr: 0.0487所以得到方程为：0082.02851.00084.00002.023++-xxx结论：在23题中计算的结果误差为4.5769，而在24中计算的结果误差为0.0487，说明对于曲线拟合来说，总会有误差，因为取样点并不是都过拟合的曲线的。
数值分析2024上机实验报告
数值分析2024上机实验报告数值分析是计算数学的一个重要分支，它研究如何用数值方法来解决数学问题。
在数值分析的学习过程中，学生需要通过上机实验来巩固理论知识，并学会使用相应的数值方法来解决实际问题。
本篇报告将详细介绍2024年度数值分析上机实验的内容和结果。
一、实验内容2024年度数值分析上机实验分为四个部分，分别是：方程求根、插值与拟合、数值积分和常微分方程的数值解。
1.方程求根这部分实验要求使用数值方法求解给定的非线性方程的根。
常见的数值方法有二分法、牛顿法、割线法等。
在实验过程中，我们需要熟悉这些数值方法的原理和实现步骤，并对不同方法的收敛性进行分析和比较。
2.插值与拟合这部分实验要求使用插值和拟合方法对给定的一组数据进行拟合。
插值方法包括拉格朗日插值、牛顿插值等；拟合方法包括最小二乘拟合、多项式拟合等。
在实验中，我们需要熟悉插值和拟合方法的原理和实现步骤，并对不同方法的精度和稳定性进行比较。
3.数值积分这部分实验要求使用数值方法计算给定函数的积分。
常见的数值积分方法有梯形法则、辛普森法则、龙贝格积分等。
在实验过程中，我们需要熟悉这些数值积分方法的原理和实现步骤，并对不同方法的精度和效率进行比较。
4.常微分方程的数值解这部分实验要求使用数值方法求解给定的常微分方程初值问题。
常见的数值方法有欧拉法、改进的欧拉法、四阶龙格-库塔法等。
在实验中，我们需要熟悉这些数值解方法的原理和实现步骤，并对不同方法的精度和稳定性进行比较。
二、实验结果在完成2024年度数值分析上机实验后，我们得到了以下实验结果：1.方程求根我们实现了二分法、牛顿法和割线法，并对比了它们的收敛速度和稳定性。
结果表明，割线法的收敛速度最快，但在一些情况下可能会出现振荡；二分法和牛顿法的收敛速度相对较慢，但稳定性较好。
2.插值与拟合我们实现了拉格朗日插值和最小二乘拟合，并对比了它们的拟合效果和精度。
结果表明，拉格朗日插值在小区间上拟合效果较好，但在大区间上可能出现振荡；最小二乘拟合在整体上拟合效果较好，但可能出现过拟合。
数值计算方法上机实验报告
数值计算方法上机实验报告
一、实验目的
本次实验的主要目的是熟悉和掌握数值计算方法，学习梯度下降法的
原理和实际应用，熟悉Python语言的编程基础知识，掌握Python语言的
基本语法。
二、设计思路
本次实验主要使用的python语言，利用python下的numpy，matplotlib这两个工具，来实现数值计算和可视化的任务。
1. 首先了解numpy的基本使用方法，学习numpy的矩阵操作，以及numpy提供的常见算法，如矩阵分解、特征值分解等。
2. 在了解numpy的基本操作后，可以学习matplotlib库中的可视化
技术，掌握如何将生成的数据以图表的形式展示出来。
3. 接下来就是要学习梯度下降法，首先了解梯度下降法的主要原理，以及具体的实际应用，用python实现梯度下降法给出的算法框架，最终
可以达到所期望的优化结果。
三、实验步骤
1. 熟悉Python语言的基本语法。
首先是熟悉Python语言的基本语法，学习如何使用Python实现变量
定义，控制语句，函数定义，类使用，以及面向对象编程的基本概念。
2. 学习numpy库的使用方法。
其次是学习numpy库的使用方法，学习如何使用numpy库构建矩阵，学习numpy库的向量，矩阵操作，以及numpy库提供的常见算法，如矩阵分解，特征值分解等。
3. 学习matplotlib库的使用方法。
数值分析第一次上机练习实验报告
数值分析第一次上机练习实验报告一、实验目的本次实验旨在通过上机练习，加深对数值分析方法的理解，并掌握实际应用中的数值计算方法。
二、实验内容1. 数值计算的基本概念和方法在本次实验中，我们首先回顾了数值计算的基本概念和方法。
数值计算是一种通过计算机进行数值近似的方法，其包括近似解的计算、误差分析和稳定性分析等内容。
2. 方程求解的数值方法接下来，我们学习了方程求解的数值方法。
方程求解是数值分析中非常重要的一部分，其目的是找到方程的实数或复数解。
我们学习了二分法、牛顿法和割线法等常用的数值求解方法，并对它们的原理和步骤进行了理论学习。
3. 插值和拟合插值和拟合是数值分析中常用的数值逼近方法。
在本次实验中，我们学习了插值和拟合的基本原理，并介绍了常见的插值方法，如拉格朗日插值和牛顿插值。
我们还学习了最小二乘拟合方法，如线性拟合和多项式拟合方法。
4. 数值积分和数值微分数值积分和数值微分是数值分析中的两个重要内容。
在本次实验中，我们学习了数值积分和数值微分的基本原理，并介绍了常用的数值积分方法，如梯形法和辛卜生公式。
我们还学习了数值微分的数值方法，如差商法和牛顿插值法。
5. 常微分方程的数值解法常微分方程是物理和工程问题中常见的数学模型，在本次实验中，我们学习了常微分方程的数值解法，包括欧拉法和四阶龙格-库塔法。
我们学习了这些方法的步骤和原理，并通过具体的实例进行了演示。
三、实验结果及分析通过本次实验，我们深入理解了数值分析的基本原理和方法。
我们通过实际操作，掌握了方程求解、插值和拟合、数值积分和数值微分以及常微分方程的数值解法等数值计算方法。
实验结果表明，在使用数值计算方法时，我们要注意误差的控制和结果的稳定性。
根据实验结果，我们可以对计算结果进行误差分析，并选择适当的数值方法和参数来提高计算的精度和稳定性。
此外，在实际应用中，我们还需要根据具体问题的特点和条件选择合适的数值方法和算法。
四、实验总结通过本次实验，我们对数值分析的基本原理和方法有了更加深入的了解。
数值分析上机实验报告
数值分析上机实验报告摘要：本报告是对数值分析课程上机实验的总结和分析，涵盖了多种算法和数据处理方法，通过对实验结果的分析，探究了数值计算的一般过程和计算的稳定性。
1. 引言数值计算是数学的一个重要分支，广泛应用于物理、金融、工程等领域。
本次实验是对数值分析课程知识的实际应用，通过上机实现算法，探究数值计算的可靠性和误差分析。
2. 实验方法本次实验中，我们实现了多种算法，包括：（1）牛顿迭代法求方程的根；（2）高斯消元法求线性方程组的解；（3）最小二乘法拟合数据点；（4）拉格朗日插值法估计函数值；（5）梯形公式和辛普森公式求积分近似值。
对于每个算法，我们都进行了多组数值和不同参数的实验，并记录了相关数据和误差。
在实验过程中，我们着重考虑了算法的可靠性和计算的稳定性。
3. 实验结果与分析在实验中，我们得到了大量的实验数据和误差分析，通过对数据的展示和分析，我们得到了以下结论：（1）牛顿迭代法求解非线性方程的根能够对算法的初始值和迭代次数进行适当的调整，从而达到更高的稳定性和可靠性。
（2）高斯消元法求解线性方程组的解需要注意到矩阵的奇异性和精度的影响，从而保证计算的准确性。
（3）最小二乘法拟合数据点需要考虑到拟合的函数形式和数据的误差范围，采取适当的数据预处理和拟合函数的选择能够提高计算的准确性。
（4）拉格朗日插值法估计函数值需要考虑到插值点的选择和插值函数的阶数，防止出现龙格现象和插值误差过大的情况。
（5）梯形公式和辛普森公式求积分近似值需要考虑到采样密度和拟合函数的选择，从而保证计算的稳定性和收敛速度。
4. 结论通过本次实验的分析和总结，我们得到了深入的认识和理解数值计算的一般过程和算法的稳定性和可靠性，对于以后的数值计算应用也提供了一定的指导和参考。
数值代数上机实验报告
试验项目名称：平方根法与改进平方根法实验内容：先用你熟悉的计算机语言将平方根法和改进平方根法编写成通用的子程序，然后用你编写的程序求解对称正定方程组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.000012.0000-69.000068.000021.0000-50.7188-8.7500-8.0000 112.00006.0000 -68.750022.000044.0000 -28.00008.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.3285-0.44380.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.00120.0139-0.09540.4208-1.21012.0624-1.0394-3.33436.2567-0.2463-7.45942.80303.69900.7277-1.7484-0.4854-3.60100.25325.18621.44100.8738-4.56541.04224.0920-2.7764-2.2148-0.89530.36654.89671.04160.1281-4.3387-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);v=v/v(1);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<mA(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)<norm(Newtonfun(X1),2) %判断先后两次迭代的大小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。
数值代数(第08周)实验报告()
A=sparse(i,j,s,n,n,3*n-2);
spy(A);
x0=A\b;
[x]=Zuiganfa(a,c,d,b);
ifnorm(x'-x0) > 1.0e-6
error('误差太大');
end
2.2.3实验（或测试）结果
得分：
??? Error using ==> Untitled2 at 14
(2)计算矩阵 的克拉特(Crout)分解
其中
【注：在实际计算中，将数组alpha和beta分别保存在数组b和c中.】
(3)解
(4)解
2.2程序清单
2.2.1写出追赶法求解三对角方程组的函数代码
得分：
function[x]=Zuiganfa(a,c,d,b);
n=length(a);
n1=length(c);
数值代数（第08周）实验报告
学号
20150
姓名
实验目的
追赶法实验
1、实验要求
给定一个 级主对角占优的三对角矩阵 ，和一个 的列向量 ，用追赶法求解线性方程组 ，其中
注：尝试用规模尽量大（即n尽量大）的矩阵来完成实验.
2、追赶法
2.1算法描述
得分：
（写出追赶法求解三对角方程组的计算过程）
(1)将矩阵 的三条对角线和常数向量 分别记为数组
%U(i,i+1)=q(i);
%end
end
2.2.2实验（或测试）代码
得分：
n=100000;
c=rand(1,n)*100;
d=rand(1,n)*100;
a = rand(1,n)*1000+c+d;
数值分析上机实验报告
数值分析上机实验报告导言：本次上机实验主要是针对数值分析课程中的一些基本算法进行实验验证。
实验内容包括迭代法、插值法、数值积分和常微分方程的数值解等。
在实验过程中，我们将会使用MATLAB进行算法的实现，并对结果进行分析。
一、迭代法迭代法是解决函数零点、方程解等问题的常用方法。
我们将选择几个常见的函数进行迭代求根的实验。
（1）二分法二分法是一种简单而有效的迭代求根法。
通过函数在区间两个端点处的函数值异号来确定函数在区间内存在零点，并通过不断缩小区间来逼近零点。
（2）牛顿法牛顿法利用函数的一阶导数和二阶导数的信息来逼近零点。
通过不断迭代更新逼近值，可以较快地求得零点。
实验结果表明，对于简单的函数，这两种迭代法都具有很好的收敛性和稳定性。
但对于一些复杂的函数，可能会出现迭代失效或者收敛速度很慢的情况。
二、插值法插值法是在给定一些离散数据点的情况下，通过构造一个插值函数来逼近未知函数的值。
本实验我们将使用拉格朗日插值和牛顿插值两种方法进行实验。
（1）拉格朗日插值拉格朗日插值通过构造一个多项式函数来逼近未知函数的值。
该多项式经过离散数据点，并且是唯一的。
该方法简单易懂，但插值点越多，多项式次数越高，插值函数的精度也就越高。
（2）牛顿插值牛顿插值利用差商的概念，通过构造一个插值多项式来逼近未知函数的值。
与拉格朗日插值相比，牛顿插值的计算过程更加高效。
但同样要求插值点的选择要合理，否则可能出现插值函数不收敛的情况。
实验结果表明，这两种插值方法都能够很好地逼近未知函数的值。
插值点的选择对插值结果有很大的影响，过多或者过少的插值点都可能导致插值结果偏离真实函数的值。
三、数值积分数值积分是一种将定积分问题转化为数值求和的方法。
本实验我们将使用复合梯形求积法和复合辛普森求积法进行实验。
（1）复合梯形求积法复合梯形求积法将定积分区间等分为若干小区间，然后使用梯形公式对每个小区间进行近似求积，最后将结果相加得到整个定积分的近似值。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
数值代数课程设计实验报告姓名： 班级： 学号： 实验日期：一、实验名称 代数的数值解法 二、实验环境 MATLAB7.0实验一、平方根法与改进平方根法一、实验要求：用熟悉的计算机语言将不选主元和列主元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:39b1(i)=15;end>> b1(40)=14;>> for i=2:83b2(i)=15;end>> b2(40)=14;>> for i=2:119b1(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 =1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.0000列主元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，二、算法描述及实验步骤：平方根法函数程序如下：1、求A 的Cholesky 分解：L L A T=; 2、求解b y =L 得y ；3、求解y x =TL 得x ；改进平方根法函数程序如下：1、求A 的Cholesky 分解：T=LDL A ; 2、求解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 =1.0000 1.00001.0000 1.00001.0000 1.00001.0000 1.00001.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.0000改进平方根法解得的解即为x2 =1.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.0000 1.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00001.00000.99991.00090.99081.09080.1010四、源程序：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 -0.75 -0.5 0 0.25 0.5 0.75iy 1.00 0.8125 0.75 1.00 1.3125 1.75 2.3125i二、算法描述及实验步骤：QR分解求解程序如下：1、求A 的QR 分解;2、计算b c 11T=Q ;3、求解上三角方程1c x =R 得x ；三、调试过程及实验结果：>> t=[-1 -0.75 -0.5 0 0.25 0.5 0.75];>> y=[ 1.00 0.8125 0.75 1.00 1.3125 1.75 2.3125]; >> 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.75 -0.5 0 0.25 0.5 0.75]; >> 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=[ 1.00 0.8125 0.75 1.00 1.3125 1.75 2.3125]; >> y=[ 1.00 0.8125 0.75 1.00 1.3125 1.75 2.3125]; >> 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 =0.3081 0.8587 1.4018 f =924/2999*x^2+10301/11996*x+4204/2999 故所求的拟合曲线为2()0.30810.8581 1.4018f x x x =++四、源程序：>> t=[-1 -0.75 -0.5 0 0.25 0.5 0.75];>> y=[ 1.00 0.8125 0.75 1.00 1.3125 1.75 2.3125]; >> plot(t,y,'r*');>> legend('实验数据(ti,yi)'); >> xlabel('t'), ylabel('y');>> title('二次多项式拟合的数据点(ti,yi)的散点图'); >> syms a b c>> t=[-1 -0.75 -0.5 0 0.25 0.5 0.75]; >> 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=[ 1.00 0.8125 0.75 1.00 1.3125 1.75 2.3125]; >> y=[ 1.00 0.8125 0.75 1.00 1.3125 1.75 2.3125]; >> 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 =0.3081 0.8587 1.4018f =924/2999*x^2+10301/11996*x+4204/2999>>word文档可自由复制编辑。