MATLAB数值分析实验三(线性方程求解及精度分析)
数值分析中求解线性方程组的MATLAB程序(6种)
![数值分析中求解线性方程组的MATLAB程序(6种)](https://img.taocdn.com/s3/m/9c7943ea81c758f5f61f67a9.png)
数值分析中求解线性方程组的MATLAB程序(6种)1.回溯法(系数矩阵为上三角)function X=uptrbk(A,B)%求解方程组,首先化为上三角,再调用函数求解[N,N]=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.'break;endfor 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);endendD=Aug;X=backsub(Aug(1:N,1:N),Aug(1:N,N+1));2.系数矩阵为下三角function x=matrix_down(A,b)%求解系数矩阵是下三角的方程组n=length(b);x=zeros(n,1);x(1)=b(1)/A(1,1);for k=2:1:nx(k)=(b(k)-A(k,1:k-1)*x(1:k-1))/A(k,k);end3.普通系数矩阵(先化为上三角,在用回溯法)function X=uptrbk(A,B)%求解方程组,首先化为上三角,再调用函数求解[N,N]=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.'break;endfor 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);endendD=Aug;X=backsub(Aug(1:N,1:N),Aug(1:N,N+1));4.三角分解法function [X,L,U]=LU_matrix(A,B)%A是非奇异矩阵%AX=B化为LUX=B,L为下三角,U为上三角%程序中并没有真正解出L和U,全部存放在A中[N,N]=size(A);X=zeros(N,1);Y=zeros(N,1);C=zeros(1,N);R=1:N;for p=1:N-1[max1,j]=max(abs(A(p:N,p)));C=A(p,:);A(p,:)=A(j+p-1,:);A(j+p-1,:)=C;d=R(p);R(p)=R(j+p-1);R(j+p-1)=d;if A(p,p)==0'A is singular.No unique solution'break;endfor k=p+1:Nmult=A(k,p)/A(p,p);A(k,p)=mult;A(k,p+1:N)=A(k,p+1:N)-mult*A(p,p+1:N);endendY(1)=B(R(1));for k=2:NY(k)=B(R(k))-A(k,1:k-1)*Y(1:k-1);endX(N)=Y(N)/A(N,N);for k=N-1:-1:1X(k)=(Y(k)-A(k,k+1:N)*X(k+1:N))/A(k,k);endL=tril(A,-1)+eye(N)U=triu(A)5.雅克比迭代法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';k6.盖斯迭代法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)break;endendX=X';k。
MATLAB计算方法3解线性方程组计算解法
![MATLAB计算方法3解线性方程组计算解法](https://img.taocdn.com/s3/m/1e3dc8c270fe910ef12d2af90242a8956becaac7.png)
MATLAB计算方法3解线性方程组计算解法线性方程组是数学中的一个重要问题,解线性方程组是计算数学中的一个基本计算,有着广泛的应用。
MATLAB是一种功能强大的数学软件,提供了多种解线性方程组的计算方法。
本文将介绍MATLAB中的三种解线性方程组的计算方法。
第一种方法是用MATLAB函数“linsolve”解线性方程组。
该函数使用高斯消元法和LU分解法求解线性方程组,可以处理单个方程组以及多个方程组的情况。
使用该函数的语法如下:X = linsolve(A, B)其中A是系数矩阵,B是常数向量,X是解向量。
该函数会根据A的形式自动选择求解方法,返回解向量X。
下面是一个使用“linsolve”函数解线性方程组的例子:A=[12;34];B=[5;6];X = linsolve(A, B);上述代码中,A是一个2×2的系数矩阵,B是一个2×1的常数向量,X是一个2×1的解向量。
运行代码后,X的值为[-4.0000;4.5000]。
第二种方法是用MATLAB函数“inv”求解逆矩阵来解线性方程组。
当系数矩阵A非奇异(可逆)时,可以使用逆矩阵求解线性方程组。
使用“inv”函数的语法如下:X = inv(A) * B其中A是系数矩阵,B是常数向量,X是解向量。
该方法先计算A的逆矩阵,然后将逆矩阵与B相乘得到解向量X。
下面是一个使用“inv”函数解线性方程组的例子:A=[12;34];B=[5;6];X = inv(A) * B;上述代码中,A是一个2×2的系数矩阵,B是一个2×1的常数向量,X是一个2×1的解向量。
运行代码后,X的值为[-4.0000;4.5000]。
第三种方法是用MATLAB函数“mldivide”(或“\”)求解线性方程组。
该函数使用最小二乘法求解非方阵的线性方程组。
使用“mldivide”函数的语法如下:X=A\B其中A是系数矩阵,B是常数向量,X是解向量。
数值分析 第三章 基于MATLAB的科学计算—线性方程组1概述
![数值分析 第三章 基于MATLAB的科学计算—线性方程组1概述](https://img.taocdn.com/s3/m/bab73266b307e87101f696bf.png)
科学计算—理论、方法及其基于MATLAB的程序实现与分析 三、 解线性方程组(线性矩阵方程)解线性方程组是科学计算中最常见的问题。
所说的“最常见”有两方面的含义:1) 问题的本身是求解线性方程组;2) 许多问题的求解需要或归结为线性方程组的求解。
关于线性方程组B A x B Ax 1-=⇒=(1)其求解方法有两类:1) 直接法:高斯消去法(Gaussian Elimination ); 2) 间接法:各种迭代法(Iteration )。
1、高斯消去法1) 引例考虑如下(梯形)线性方程组:()⎪⎩⎪⎨⎧==+==+-=⇒⎪⎪⎪⎭⎫ ⎝⎛=⎪⎪⎪⎭⎫ ⎝⎛⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡--⇔⎪⎩⎪⎨⎧==-=+-5.0141315.3221122004301211214322332321321332321x x x x x x x x x x x x x x x 高斯消去法的求解思路:把一般的线性方程组(1)化成(上或下)梯形的形式。
2)高斯消去法——示例考虑如下线性方程组:⎪⎪⎪⎭⎫⎝⎛-=⎪⎪⎪⎭⎫ ⎝⎛⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡---⇔⎪⎩⎪⎨⎧=++-=-+-=+-306015129101.2001.221113060129501.2001.221321321321321x x x x x x x x x x x x1) 第一个方程的两端乘12加到第二个方程的两端,第一个方程的两端乘-1加到第三个方程的两端,得⎪⎪⎪⎭⎫⎝⎛-=⎪⎪⎪⎭⎫ ⎝⎛⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡--3060031110001.0001.00111321x x x2) 第二个方程的两端乘001.010-加到第三个方程的两端,得 ⎪⎪⎪⎭⎫⎝⎛-=⎪⎪⎪⎭⎫ ⎝⎛⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡--60600311010001.0001.00111321x x x3) 从上述方程组的第三个方程依此求解,得()⎪⎩⎪⎨⎧==+-==+-=600300001.03100024011332321x x x x x x 3)高斯消去法的不足及其改进——高斯(全、列)主元素消去法在上例中,由于建模、计算等原因,系数2.001而产生0.0005的误差,实际求解的方程组为⎪⎪⎪⎭⎫⎝⎛-=⎪⎪⎪⎭⎫ ⎝⎛⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡---306015129101.20005.22111321x x x ⎪⎩⎪⎨⎧===⇒70.4509.30142.2565321x x x 注:数值稳定的算法高斯列主元素消去法就是在消元的每一步选取(列)主元素—一列中绝对值最大的元取做主元素,高斯列主元素消去法是数值稳定的方法。
matlab线性方程组数值求解实验报告
![matlab线性方程组数值求解实验报告](https://img.taocdn.com/s3/m/2d9967c6240c844768eaee17.png)
湖南大学电气与信息工程学院 《数值计算》课程 上机实验报告一. 实验目的:了解gauss 消去法和迭代法matlab 算法实现求任意方程组的根。
二. 实验内容:用gauss 消去法和迭代法求解下列线性方程组:263234323923321321321=++=++=++x x x x x x x x x1.求出gauss 消去法的上三角矩阵和方程组的解321,,x x x ,并在命令窗口显示;2.显示迭代法求解过程中所有结果(,,,,,,,,,321131*********NN N x x x x x x x x x ⋯⋯)要求求解精度达到10^-5.三. 算法介绍或方法基础1) 消去法:消元过程:设0)0(11≠a ,令乘数)0(11)0(11/a a m i i -=,做(消去第i 个方程组的i x )操作1i m ×第1个方程+第i 个方程(i=2,3,.....n )则第i 个方程变为1)1(2)1(2...i n in i b x a x a =++ 这样消去第2,3,。
,n 个方程的变元i x 后。
原线性方程组变为:⎪⎪⎪⎩⎪⎪⎪⎨⎧=++=++=++)1()1(2)1(2)1(2)1(22)1(22)0(1)0(11)0(11... . .... ...n n nn n n n n n b x a x a b x a x a b x a x a 这样就完成了第1步消元。
回代过程:在最后的一方程中解出n x ,得:)1()1(/--=n nn n n n a b x再将n x 的值代入倒数第二个方程,解出1-n x ,依次往上反推,即可求出方程组的解: 其通项为3,...1-n 2,-n k /)()1(1)1()1(=-=-+=--∑k kk nk j j k kj k kk a x abx高斯赛德尔迭代法:由雅可比迭代公式可知,在迭代的每一步计算过程中是用()k x的全部分量来计算()1+k x的所有分量,显然在计算第i 个分量()1+k ix 时,已经计算出的最新分量()()1111+-+k i k x ,...,x 没有被利用,从直观上看,最新计算出的分量可能比旧的分量要好些.因此,对这些最新计算出来的第1+k 次近似()1+k x 的分量()1+k jx 加以利用,就得到所谓解方程组的高斯—塞德(Gauss-Seidel )迭代法.把矩阵A 分解成U L D A --= (6)其中()nn a ,...,a ,a diag D 2211=,U ,L --分别为A 的主对角元除外的下三角和上三角部分,于是,方程组(1)便可以写成 ()b Ux x L D +=-即 22f x B x +=其中()()b L D f ,U L D B 1212---=-= (7)以2B 为迭代矩阵构成的迭代法(公式)()()221f x B x k k +=+ (8)称为高斯—塞德尔迭代法(公式),用 量表示的形式为⎩⎨⎧[],...,,k ,n ,,i x a x a b a xi j n i j )k (j ij )k (j ij i ii)k (i21021111111==∑∑--=-=+=++ (9)由此看出,高斯—塞德尔迭代法的一个明显的优点是,在电算时,只需一组存储单元(计算出()1+k ix 后()k ix 不再使用,所以用()1+k ix 冲掉()k ix ,以便存放近似解.四.程序1)消去法:function x=gauss(A,b)n=length(b);A=[A,b];for k=1:(n-1)A((k+1):n,(k+1):(n+1))=A((k+1):n,(k+1):(n+1))... -A((k+1):n,k)/A(k,k)*A(k,(k+1):(n+1));A((k+1):n,k)=zeros(n-k,1);Aendx=zeros(n,1);x(n)=A(n,n+1)/A(n,n);for k=n-1:-1:1x(k,:)=(A(k,n+1)-A(k,(k+1):n)*x((k+1):n))/A(k,k);end2)迭代法:function EX()a=input('请输入系数矩阵a:');b=input('请输入矩阵b:');N=input('请输入最大迭代次数N:');esp=input('请输入近似解的误差限:');if any(diag(a))==0error('系数矩阵错误,迭代终止!')endD=diag(diag(a));X0=zeros(size(b));x1=0;x2=0;x3=0;X1=[x1;x2;x3];h=inv(D)*b;B=inv(D)*(D-a);B1=triu(B);B2=tril(B);k=1;fprintf('高斯-赛德尔迭代法');fprintf('第0次迭代得:')disp(X1');while k<=Nx1=h(1,1)+B1(1,:)*X0;X1=[x1;x2;x3];x2=h(2,1)+B1(2,:)*X0+B2(2,:)*X1;X1=[x1;x2;x3];x3=h(3,1)+B2(3,:)*X1;X1=[x1;x2;x3];if norm(X1-X0,inf)<espfprintf('已满足误差限。
MATLAB求解方程解析解和数值解
![MATLAB求解方程解析解和数值解](https://img.taocdn.com/s3/m/d64e6d1b52ea551810a68760.png)
辽宁工程技术大学上机实验报告用MATLAB求解质点振动方程振动是日常生活和工程技术中常见的一种运动形式。
利用常系数线性微分方程的理论来讨论有关自由振动和强迫振动的相关问题。
利用MA TLAB数学软件大致可分四类情况:(1)无阻尼自由振动情况;(2)有阻尼自由振动;(3)无阻尼强迫振动;(4)有阻尼强迫振动求其数值解和解析解;MATLAB软件求解微分方程解析解的命令“dsolve()”求通解的命令格式:(’微分方程’,’自变量’)注:微分方程在输入时,一阶导数y’应输入Dy,y’’应输入D2y等,D应大写。
1,无阻尼自由振动情况:常见的数学摆的无阻尼微小振动方程代码如下:>> t=0:pi/50:2*pi;>> y=2*sin(3*t+2);>>plot(t,y,'b')2,有阻尼自由振动由无阻尼振动的通解可以看出,无阻尼振动是按照正弦规律运动的,摆动似乎可以无限期的进行下去,但事实上,空气从在阻力,在运动时,我们必须把空气阻力考虑在内,所以我们得到有阻尼摆动方程为:记u/m=2n,g/l=w^2,这里n,w是正常数,所以:y=dsolve('D2y+2*n*Dy+w^2*y=0','t'); (4.43)解得:y = C3*exp(-t*(n + ((n + w)*(n - w))^(1/2))) + C2*exp(-t*(n - ((n + w)*(n - w))^(1/2)))(1)小阻尼情形:n<w时,方程(4.43)的通解为:y=exp(-n*t)*(c1*cos(w1*t)+c2*sin(w1*t))和前面无阻尼的情形一样,可以把上式的通解改写为一下形式:y=A*exp(-n*t)*sin(w1*t+Q), (4.45)这里的A,Q为任意常数。
用matlab 操作得到:t=0:0.1:10;y=3*exp(-0.1*t).*sin(5*t+4);plot(t,y,'k-')如图:由(4.45)可见,摆动的运动不是周期的,振动的幅度随着时间的增加而不断减小。
MATLAB方程组求解实验报告
![MATLAB方程组求解实验报告](https://img.taocdn.com/s3/m/fd90a35233687e21af45a9f1.png)
八、教师评语
签名:
日期: 年 月 日
成绩
14 -2 2
0 33 8
-2 7 38
3)
a=[2 1 1;1 3 1;1 1 4];
b=[12 -3 1;-1 30 7;-3 6 34];
>> a-b
ans =
-10 4 0
2 -27 -6
4 -5 -30
4)
a=[2 1 1;1 3 1;1 1 4];
b=[12 -3 1;-1 30 7;-3 6 34];
方程A*X=b变形成QRX=b则X=R\(Q\b)
命令[Q,R]=qr(A)
2求线性齐次方程组的通解
在Matlab中,函数null用来求解零空间,即满足A·X=0的解空间,实际上是求出解空间的一组基(基础解系)。
格式:z=null %z的列向量为方程组的正交规范基,满足
% z的列向量是方程AX=0的有理基
1)
a=[2 1 1;1 3 1;1 1 4];
b=[12 -3 1;-1 30 7;-3 6 34];
>> a*b
ans =
20 30 43
6 93 56
-1 51 144
2)
a=[2 1 1;1 3 1;1 1 4];
b=[12 -3 1;-1 30 7;-3 6 34];
>> a+b
ans =
-2.3529 12.5294 0.7059
-2.2353 0.3529 9.4706
数值分析matlab实验报告
![数值分析matlab实验报告](https://img.taocdn.com/s3/m/284fa9c303d276a20029bd64783e0912a2167cee.png)
数值分析matlab实验报告《数值分析MATLAB实验报告》摘要:本实验报告基于MATLAB软件进行了数值分析实验,通过对不同数学问题的数值计算和分析,验证了数值分析方法的有效性和准确性。
实验结果表明,MATLAB在数值分析领域具有较高的应用价值和实用性。
一、引言数值分析是一门研究利用计算机进行数值计算和分析的学科,其应用范围涵盖了数学、物理、工程等多个领域。
MATLAB是一种常用的数值计算软件,具有强大的数值分析功能,能够进行高效、准确的数值计算和分析,因此在科学研究和工程实践中得到了广泛的应用。
二、实验目的本实验旨在通过MATLAB软件对数值分析方法进行实验验证,探究其在不同数学问题上的应用效果和准确性,为数值分析方法的实际应用提供参考和指导。
三、实验内容1. 利用MATLAB进行方程求解实验在该实验中,利用MATLAB对给定的方程进行求解,比较数值解和解析解的差异,验证数值解的准确性和可靠性。
2. 利用MATLAB进行数值积分实验通过MATLAB对给定函数进行数值积分,比较数值积分结果和解析积分结果,验证数值积分的精度和稳定性。
3. 利用MATLAB进行常微分方程数值解实验通过MATLAB对给定的常微分方程进行数值解,比较数值解和解析解的差异,验证数值解的准确性和可靠性。
四、实验结果与分析通过对以上实验内容的实际操作和分析,得出以下结论:1. 在方程求解实验中,MATLAB给出的数值解与解析解基本吻合,验证了MATLAB在方程求解方面的高准确性和可靠性。
2. 在数值积分实验中,MATLAB给出的数值积分结果与解析积分结果基本吻合,验证了MATLAB在数值积分方面的高精度和稳定性。
3. 在常微分方程数值解实验中,MATLAB给出的数值解与解析解基本吻合,验证了MATLAB在常微分方程数值解方面的高准确性和可靠性。
五、结论与展望本实验通过MATLAB软件对数值分析方法进行了实验验证,得出了数值分析方法在不同数学问题上的高准确性和可靠性。
数值分析实验报告matlab
![数值分析实验报告matlab](https://img.taocdn.com/s3/m/fada48b49f3143323968011ca300a6c30c22f1bc.png)
数值分析实验报告matlab数值分析实验报告引言:数值分析是一门研究利用计算机数值方法解决数学问题的学科,它在科学计算、工程设计、金融分析等领域具有重要的应用价值。
本实验报告旨在通过使用MATLAB软件,探索数值分析的基本原理和方法,并通过实际案例加深对数值分析的理解。
一、误差分析在数值计算中,误差是无法避免的。
误差分析是数值分析中的重要一环,它帮助我们了解数值计算的准确性和稳定性。
在实验中,我们通过计算机模拟了一个简单的数学问题,并分别计算了绝对误差和相对误差。
通过比较不同算法的误差大小,我们可以选择最适合的算法来解决实际问题。
二、插值与拟合插值和拟合是数值分析中常用的方法,它们可以通过已知的数据点来推导出未知数据点的近似值。
在本实验中,我们通过MATLAB的插值函数和拟合函数,分别进行了插值和拟合的实验。
通过比较不同插值和拟合方法的结果,我们可以选择最适合的方法来处理实际问题。
三、数值积分数值积分是数值分析中的重要内容,它可以用来计算曲线下的面积或函数的积分值。
在实验中,我们通过MATLAB的数值积分函数,对一些简单的函数进行了积分计算。
通过比较数值积分和解析积分的结果,我们可以评估数值积分的准确性和稳定性,并选择最适合的积分方法来解决实际问题。
四、常微分方程的数值解法常微分方程是数值分析中的重要内容,它可以用来描述许多自然现象和工程问题。
在实验中,我们通过MATLAB的常微分方程求解函数,对一些简单的微分方程进行了数值解法的计算。
通过比较数值解和解析解的结果,我们可以评估数值解法的准确性和稳定性,并选择最适合的数值解法来解决实际问题。
五、线性方程组的数值解法线性方程组是数值分析中的经典问题,它在科学计算和工程设计中广泛应用。
在实验中,我们通过MATLAB的线性方程组求解函数,对一些简单的线性方程组进行了数值解法的计算。
通过比较数值解和解析解的结果,我们可以评估数值解法的准确性和稳定性,并选择最适合的数值解法来解决实际问题。
MATLAB试验报告线性方程组求解
![MATLAB试验报告线性方程组求解](https://img.taocdn.com/s3/m/f030dd4759fb770bf78a6529647d27284b7337b0.png)
「-「
f3)
f0)
f2)
f-1)
x=
2
,x=
4
,x=
10
,x=
7
,x=
3
1
3
2
2
3
11
4
3
5
2
\^)
V37
R3的一组基。
生成R3,从x,x,x,x,x中找出
"-4
设A = 7
11
3 12、
-11 0.分别用poly、roots和eig计算A的特征值.
12 3,
7.判断下列矩阵是否相似(提示:用Jordan标准型)?
,
壬式。
及rref函数,计算下列线性方程组的
x-3x=0
5x+x=8b)
-x+4x=2
2 3
:1
4,%=10
27 V117
是线性相关的.
戋性无关的
生方程组的通解。
至解
=2
二10
8
来计算A-1的第三列,并将结果
通解
2x+x-3x=0
4x+5x+x=8
2x+4x+4x=8
1 1 2 3
实验序号:实验2
实验项目名称:线性方程组求解
>> x1=det(D1)/det(D);x2=det(D2)/det(D);
>> x1,x2
x1 =
-1
Байду номын сангаас>> a=[2,3,-1;5,1,4];
>> rref(a)
f2 0 0)
f2 0
A=
0 4 0
数值分析实验报告三
![数值分析实验报告三](https://img.taocdn.com/s3/m/b3283a287f21af45b307e87101f69e314332fafe.png)
grid
[k,x,wuca,yx]=erfen (﹣1,1,10^-5)
2)运行结果
ans =
0 -1.0000 1.0000 0 1.0000 -11.6321 10.7183 -1.0000
ans =
1.0000 0 1.0000 0.5000 0.5000 -1.0000 10.7183 4.6487
ans =
11.0000 0.0898 0.0908 0.0903 0.0005 -0.0076 0.0033 -0.0021
ans =
12.0000 0.0903 0.0908 0.0906 0.0002 -0.0021 0.0033 0.0006
ans =
13.0000 0.0903 0.0906 0.0905 0.0001 -0.0021 0.0006 -0.0008
ans =
7.0000 0.1256 0.0008 0.0033 0.0262
ans =
8.0000 0.1240 0.0002 0.0016 0.0129
ans =
9.0000 0.1233 0.0000 0.0007 0.0056
ans =
9.0000 0.1233 0.0000 0.0007 0.0056
(2)、Use the iteration method ,the initial value .
2、The equation has two roots near 0.1.
Determine them by means ofNewton’s method.
(with accuracy )
3、用迭代法求方程 附近的一个根。方程写成下
k = 9
matlab线性方程组求解实验报告
![matlab线性方程组求解实验报告](https://img.taocdn.com/s3/m/6196cb3087c24028915fc318.png)
4. (1)Guass 消去法 请输入系数方阵 A=[3 2 1;2 3 1;1 2 3] 请输入列向量 B=[39;34;26] 合并后的增广矩阵为 c= 3 2 1 2 3 2 1 1 3 39 34 26
开始进行列主元消元法 c=
3 2 1
2 3 2
1 1 3
39 34 26
c= 3.00000000000000 0 0 消元后的矩阵 C c= 3.00000000000000 0 0 2.00000000000000 1.66666666666667 0 1.00000000000000 39.00000000000000 0.33333333333333 8.00000000000000 2.40000000000000 6.60000000000000 2.00000000000000 1.66666666666667 1.33333333333333 1.00000000000000 39.00000000000000 0.33333333333333 8.00000000000000 2.66666666666667 13.00000000000000
A=[2/3 1/3;2/3 1/3;1/3 2/3];%第二次试验的方程组
B=[39/3;34/3;26/3]; c=[];d=[];e=[]; x=zeros(1,3);y=zeros(1,3); y(1)=B(1)-A(1,1)*x(2)-A(1,2)*x(3); y(2)=B(2)-A(2,1)*x(1)-A(2,2)*x(3); y(3)=B(3)-A(3,1)*x(1)-A(3,2)*x(2); c=[c;y(1)];d=[d;y(2)];e=[e;y(3)]; while i<20 %abs(x(1)-y(1))>0.001&abs(x(2)-y(2))>0.001&abs(x(3)-y(3))>0.001 x=y; %x y(1)=B(1)-A(1,1)*x(2)-A(1,2)*x(3);c=[c;y(1)]; y(2)=B(2)-A(2,1)*x(1)-A(2,2)*x(3);d=[d;y(2)]; y(3)=B(3)-A(3,1)*x(1)-A(3,2)*x(2);e=[e;y(3)]; i=i+1; end plot(c,'-r'); hold on; plot(d,'-k'); plot(e,'-b'); hold off legend('x1','x2','x3','Location','NorthWest');
数值分析实验报告(Matlab实现)(同名6593)
![数值分析实验报告(Matlab实现)(同名6593)](https://img.taocdn.com/s3/m/fb4e6bfa581b6bd97e19ea9d.png)
学生实验报告实验课程名称数值分析开课实验室数学与统计学院实验室学院2010 年级数学与应用数学专业班01班学生姓名学号开课时间2012 至2013 学年第一学期end y=x;format short ;% 设置为默认格式显示,显示5位(2) 建立MATLAB 界面利用MA TLAB 的GUI 建立如下界面求解线性方程组:详见程序。
五、 计算实例、数据、结果、分析下面我们对以上的结果进行测试,求解:⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡-=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡------725101391444321131243301024321x x x x 输入数据后点击和,得到如下结果:更改以上数据进行测试,求解如下方程组:123443211343212343112341x x x x ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥-⎢⎥⎢⎥⎢⎥-⎣⎦⎣⎦⎣⎦得到如下结果:六、 实验中遇到的问题及解决办法在本实验中,遇到的问题主要有两个:(1) 如何将上述的Gauss 消元法的算法在MA TLAB 中实现针对此问题我借鉴了网上以及 课本上的算法的MATLAB 实现的程序;(2) 如何将建立界面使得可以随意输入想要求解的相关矩阵后就可以直接求解针对此问题,我通过网上的一些关于MA TLAB 的GUI 设计的相关资料,总结经验完成了此项任务。
七、 实验结论通过以上的测试,我们发现以上算法和程序能够求出线性方程组的比较精确解。
八、 参考文献[1]杨大地,王开荣 .2006.数值分析.北京:科学出版社[2]何光辉.2008. 数值分析实验. 重庆大学数理学院数学实验教学中心 [3]百度文库,百度知道教师签名年 月 日详见程序。
五、计算实例、数据、结果、分析下面我们对以上的问题进行测试:输入数据:计算结果如下:当x=2.101时,x=4.234时,同理可以测试(4)中的5的值。
六、实验中遇到的问题及解决办法在本实验中,遇到的问题主要有两个:(3)如何将上述的插值的算法在MA TLAB中实现针对此问题我借鉴了网上以及课本上的算法的MATLAB实现的程序;(4)如何将建立界面使得可以随意输入想要求解的相关矩阵后就可以直接求解针对此问题,我通过网上的一些关于MA TLAB的GUI设计的相关资料,总结经验完成了此项任务。
数值分析MATLAB科学计算—线性方程组
![数值分析MATLAB科学计算—线性方程组](https://img.taocdn.com/s3/m/98cb2d5cf46527d3240ce054.png)
科学计算—理论、方法及其基于MATLAB的程序实现与分析 三、 解线性方程组(线性矩阵方程)解线性方程组是科学计算中最常见的问题。
所说的“最常见”有两方面的含义:1) 问题的本身是求解线性方程组;2) 许多问题的求解需要或归结为线性方程组的求解。
关于线性方程组B A x B Ax 1-=⇒=(1)其求解方法有两类:1) 直接法:高斯消去法(Gaussian Elimination ); 2) 间接法:各种迭代法(Iteration )。
1、高斯消去法1) 引例考虑如下(梯形)线性方程组:()⎪⎩⎪⎨⎧==+==+-=⇒⎪⎪⎪⎭⎫ ⎝⎛=⎪⎪⎪⎭⎫ ⎝⎛⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡--⇔⎪⎩⎪⎨⎧==-=+-5.0141315.3221122004301211214322332321321332321x x x x x x x x x x x x x x x 高斯消去法的求解思路:把一般的线性方程组(1)化成(上或下)梯形的形式。
2)高斯消去法——示例考虑如下线性方程组:⎪⎪⎪⎭⎫ ⎝⎛-=⎪⎪⎪⎭⎫ ⎝⎛⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡---⇔⎪⎩⎪⎨⎧=++-=-+-=+-306015129101.2001.221113060129501.2001.221321321321321x x x x x x x x x x x x 1) 第一个方程的两端乘12加到第二个方程的两端,第一个方程的两端乘-1加到第三个方程的两端,得⎪⎪⎪⎭⎫⎝⎛-=⎪⎪⎪⎭⎫ ⎝⎛⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡--3060031110001.0001.00111321x x x2) 第二个方程的两端乘001.010-加到第三个方程的两端,得 ⎪⎪⎪⎭⎫⎝⎛-=⎪⎪⎪⎭⎫ ⎝⎛⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡--60600311010001.0001.00111321x x x3) 从上述方程组的第三个方程依此求解,得()⎪⎩⎪⎨⎧==+-==+-=600300001.03100024011332321x x x x x x 3)高斯消去法的不足及其改进——高斯(全、列)主元素消去法在上例中,由于建模、计算等原因,系数2.001而产生0.0005的误差,实际求解的方程组为⎪⎪⎪⎭⎫⎝⎛-=⎪⎪⎪⎭⎫ ⎝⎛⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡---306015129101.20005.22111321x x x ⎪⎩⎪⎨⎧===⇒70.4509.30142.2565321x x x注:数值稳定的算法高斯列主元素消去法就是在消元的每一步选取(列)主元素—一列中绝对值最大的元取做主元素,高斯列主元素消去法是数值稳定的方法。
MATLAB数值分析实验三(线性方程求解及精度分析)
![MATLAB数值分析实验三(线性方程求解及精度分析)](https://img.taocdn.com/s3/m/ec26429d84868762caaed570.png)
佛山科学技术学院实 验 报 告课程名称 数值分析实验项目 数值积分专业班级 机械工程 姓 名 余红杰 学 号 2111505010 指导教师 陈剑 成 绩 日 期 月 日一、实验目的1、 掌握程序的录入和matlab 的使用和操作;2、 了解影响线性方程组解的精度的因素——方法与问题的性态。
3、 学会Matlab 提供的“\”的求解线性方程组。
二、实验要求1、按照题目要求完成实验内容;2、写出相应的Matlab 程序;3、给出实验结果(可以用表格展示实验结果);4、分析和讨论实验结果并提出可能的优化实验。
5、写出实验报告。
三、实验步骤1、用LU 分解及列主元高斯消去法解线性方程组a)⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=⎪⎪⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎪⎪⎭⎫ ⎝⎛----15900001.582012151526099999.23107104321x x x x , 输出b Ax =中系数LU A =分解的矩阵L 和U ,解向量x 和)det(A ;用列主元法的行交换次序解向量x 和求)det(A ;比较两种方法所得结果。
2、用列主高斯消元法解线性方程组b Ax =。
(1)、⎪⎪⎪⎭⎫ ⎝⎛=⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎭⎫ ⎝⎛--11134.981.4987.023.116.427.199.103.601.3321x x x(2)、⎪⎪⎪⎭⎫ ⎝⎛=⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎭⎫ ⎝⎛--11134.981.4990.023.116.427.199.103.600.3321x x x 分别输出)det(,,A b A ,解向量x ,(1)中A 的条件数。
分析比较(1)、(2)的计算结果3、线性方程组b Ax =的A 和b 分别为⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=1095791068565778710A ,⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=31332332b 则解T x ),1,1,1,1(=. 用MATLAB 内部函数求)det(A 和A 的所有特征值和2)(A cond . 若令⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=+98.99599.6989.998.585604.508.72.71.8710A A δ, 求解b x x A A =++))((δδ,输出向量x δ和2x δ,从理论结果和实际计算两方面分析线性方程组b Ax =解的相对误差22/x xδ以及A 的相对误差/A A δ的关系。
数值分析matlab实验报告
![数值分析matlab实验报告](https://img.taocdn.com/s3/m/e2aae840a9114431b90d6c85ec3a87c241288a50.png)
数值分析matlab实验报告数值分析 Matlab 实验报告一、实验目的数值分析是研究各种数学问题数值解法的学科,Matlab 则是一款功能强大的科学计算软件。
本次实验旨在通过使用 Matlab 解决一系列数值分析问题,加深对数值分析方法的理解和应用能力,掌握数值计算中的误差分析、数值逼近、数值积分与数值微分等基本概念和方法,并培养运用计算机解决实际数学问题的能力。
二、实验内容(一)误差分析在数值计算中,误差是不可避免的。
通过对给定函数进行计算,分析截断误差和舍入误差的影响。
例如,计算函数$f(x) =\sin(x)$在$x = 05$ 附近的值,比较不同精度下的结果差异。
(二)数值逼近1、多项式插值使用拉格朗日插值法和牛顿插值法对给定的数据点进行插值,得到拟合多项式,并分析其误差。
2、曲线拟合采用最小二乘法对给定的数据进行线性和非线性曲线拟合,如多项式曲线拟合和指数曲线拟合。
(三)数值积分1、牛顿柯特斯公式实现梯形公式、辛普森公式和柯特斯公式,计算给定函数在特定区间上的积分值,并分析误差。
2、高斯求积公式使用高斯勒让德求积公式计算积分,比较其精度与牛顿柯特斯公式的差异。
(四)数值微分利用差商公式计算函数的数值导数,分析步长对结果的影响,探讨如何选择合适的步长以提高精度。
三、实验步骤(一)误差分析1、定义函数`compute_sin_error` 来计算不同精度下的正弦函数值和误差。
```matlabfunction value, error = compute_sin_error(x, precision)true_value = sin(x);computed_value = vpa(sin(x), precision);error = abs(true_value computed_value);end```2、在主程序中调用该函数,分别设置不同的精度进行计算和分析。
(二)数值逼近1、拉格朗日插值法```matlabfunction L = lagrange_interpolation(x, y, xi)n = length(x);L = 0;for i = 1:nli = 1;for j = 1:nif j ~= ili = li (xi x(j))/(x(i) x(j));endendL = L + y(i) li;endend```2、牛顿插值法```matlabfunction N = newton_interpolation(x, y, xi)n = length(x);%计算差商表D = zeros(n, n);D(:, 1) = y';for j = 2:nfor i = j:nD(i, j) =(D(i, j 1) D(i 1, j 1))/(x(i) x(i j + 1));endend%计算插值结果N = D(1, 1);term = 1;for i = 2:nterm = term (xi x(i 1));N = N + D(i, i) term;endend```3、曲线拟合```matlab%线性最小二乘拟合p = polyfit(x, y, 1);y_fit_linear = polyval(p, x);%多项式曲线拟合p = polyfit(x, y, n);% n 为多项式的次数y_fit_poly = polyval(p, x);%指数曲线拟合p = fit(x, y, 'exp1');y_fit_exp = p(x);```(三)数值积分1、梯形公式```matlabfunction T = trapezoidal_rule(f, a, b, n)h =(b a) / n;x = a:h:b;y = f(x);T = h ((y(1) + y(end))/ 2 + sum(y(2:end 1)));end```2、辛普森公式```matlabfunction S = simpson_rule(f, a, b, n)if mod(n, 2) ~= 0error('n 必须为偶数');endh =(b a) / n;x = a:h:b;y = f(x);S = h / 3 (y(1) + 4 sum(y(2:2:end 1))+ 2 sum(y(3:2:end 2))+ y(end));end```3、柯特斯公式```matlabfunction C = cotes_rule(f, a, b, n)h =(b a) / n;x = a:h:b;y = f(x);w = 7, 32, 12, 32, 7 / 90;C = h sum(w y);end```4、高斯勒让德求积公式```matlabfunction G = gauss_legendre_integration(f, a, b)x, w = gauss_legendre(5);%选择适当的节点数t =(b a) / 2 x +(a + b) / 2;G =(b a) / 2 sum(w f(t));end```(四)数值微分```matlabfunction dydx = numerical_derivative(f, x, h)dydx =(f(x + h) f(x h))/(2 h);end```四、实验结果与分析(一)误差分析通过不同精度的计算,发现随着精度的提高,误差逐渐减小,但计算时间也相应增加。
MATLAB线性方程组求解方法
![MATLAB线性方程组求解方法](https://img.taocdn.com/s3/m/40df892d5901020207409c29.png)
——GW318 物联网实验室学术活动
MATLAB 线性方程组求解方法
1.线性方程组的问题 在工程计算中,一个重要的问题是线性方程组的求解。在矩阵表示法中,上述问题可以 表述为给定两个矩阵 A 和 B,是否存在惟一的解 X 使得 AX=B 或 XA=B。 例如,求解方程 3x=6 就可以将矩阵 A 和 B 看成是标量的一种情况,最后得出该方 程的 解为 x=6/3=2。 尽管在标准的数学中并没有矩阵除法的概念,但 MATLAB 7.0 采用了与解标量方程中 类 似的约定,用除号来表示求解线性方程的解。MATLAB 7.0 采用第 2 章介绍过的运算 符斜杠 ’/’和运算符反斜杠’\ ’来表示求线性方程的解,其具体含义如下: • X=A\B 表示求矩阵方程 AX=B 的解; • X=B/A 表示求矩阵方程 XA=B 的解。 对于 X=A\B,要求矩阵 A 和矩阵 B 有相同的行数,X 和 B 有相同的列数,X 的行 数等于 矩阵 A 的列数。X=B\A 行和列的性质则与之相反。 在实际情况中,形式 AX=B 的线性方程组比形式 XA=B 的线性方程组要常见得多。 因此 反斜杠’\’用得更多。本小节的内容也主要针对反斜杠’\’除法进行介绍。斜杠’/’ 除法的性质可以 由恒等变换式得到(B/A)’=(A’\B’)。 系数矩阵 A 不一定要求是方阵,矩阵 A 可以是 m×n 的矩阵,有如下 3 种情况: • m=n 恰定方程组,MATLAB7.0 会寻求精确解; • m>n 超定方程组,MATLAB7.0 会寻求最小二乘解; • m<n 欠定方程组,MATLAB7.0 会寻求基本解,该解最多有 m 个非零元素。 值得注意的是用 MATLAB 7.0 求解这种问题时,并不采用计算矩阵的逆的方法。针对 不 同的情况,MATLAB7.0 会采用不同的算法来解线性方程组。 2.线性方程组的一般解 线性方程组 AX=B 的一般解给出了满足 AX=B 的所有解。线性方程组的一般解可以通 过 下面的步骤得到。 • 解相应的齐次方程组 AX=0, 求得基础解。 可以使用函数 null()来得到基础解。 语 句 null(A) 返回齐次方程组 AX=0 的一个基础解,其他基础解与 null(A)是线性关系。 • 求非齐次线性方程组 AX=B,得到一个特殊解。 • 非齐次线性方程组 AX=B 的一般解等于基础解的线性组合加上特殊解。 在后面的章节将介绍求非齐次线性方程组 AX=B 特殊解的方法。 3.恰定方程组的求解
MATLAB方程组求解实验报告
![MATLAB方程组求解实验报告](https://img.taocdn.com/s3/m/b9e3d823ae1ffc4ffe4733687e21af45b207fe75.png)
MATLAB方程组求解实验报告实验报告:MATLAB方程组求解一、引言在工程和科学领域的研究中,常常需要求解一系列的线性或非线性方程组。
MATLAB是一种强大的数学软件,可以用于解决方程组求解的问题。
本实验旨在通过实例介绍MATLAB求解方程组的方法和应用。
二、方程组的定义与求解方法方程组是一组包含多个未知数的方程的集合。
求解方程组即求解方程组中的未知数,使得方程组中的每个方程都成立。
对于线性方程组,可以使用矩阵表示。
例如:Ax=b其中A是一个已知的mxn的矩阵,x是待求解的向量,b是已知的向量。
MATLAB提供了多种求解线性方程组的方法,如高斯消元法、LU分解法和迭代法等。
对于非线性方程组,一般使用数值方法求解。
常见的数值方法有牛顿法、割线法和迭代法等。
MATLAB中的fzero函数可以用于求解非线性方程组。
三、MATLAB求解线性方程组的实例1.高斯消元法考虑以下线性方程组:3x+2y-z=12x-2y+4z=-2-x+0.5y-z=0可以通过高斯消元法求解该方程组。
在MATLAB中,可以使用linsolve函数进行求解。
2.LU分解法LU分解是一种常用的求解线性方程组的方法。
通过将系数矩阵分解为上三角矩阵U和下三角矩阵L的乘积来求解方程组。
在MATLAB中,可以使用lu函数进行LU分解。
四、MATLAB求解非线性方程组的实例1.牛顿法考虑以下非线性方程组:x^2+y^2=1(x-1)^2+y^2=1可以通过牛顿法求解该方程组。
在MATLAB中,可以使用fsolve函数进行求解。
2.迭代法考虑以下非线性方程组:x^2+y^2=2x+y=1可以通过迭代法求解该方程组。
在MATLAB中,设置初始值,并使用循环迭代的方法逐步逼近方程的解。
五、实验步骤和结果1.线性方程组求解构建线性方程组Ax = b,并使用linsolve函数进行求解。
2.非线性方程组求解构建非线性方程组,并使用fsolve函数进行求解。
matlab数值分析实验报告
![matlab数值分析实验报告](https://img.taocdn.com/s3/m/f7540277effdc8d376eeaeaad1f34693daef10e2.png)
matlab数值分析实验报告Matlab数值分析实验报告引言数值分析是一门研究利用计算机进行数值计算和模拟的学科,它在科学计算、工程技术和金融等领域有着广泛的应用。
本次实验报告将介绍在Matlab环境下进行的数值分析实验,包括数值微分、数值积分和线性方程组求解等内容。
一、数值微分数值微分是通过数值方法计算函数的导数,常用的数值微分方法有前向差分、后向差分和中心差分。
在Matlab中,可以使用diff函数来计算函数的导数。
例如,对于函数f(x)=x^2,在Matlab中可以使用如下代码进行数值微分的计算:```matlabsyms x;f = x^2;df = diff(f, x);```二、数值积分数值积分是通过数值方法计算函数的定积分,常用的数值积分方法有梯形法则、辛普森法则和龙贝格积分法。
在Matlab中,可以使用trapz、quad和integral等函数来进行数值积分的计算。
例如,对于函数f(x)=sin(x),可以使用如下代码进行数值积分的计算:```matlabx = linspace(0, pi, 100);y = sin(x);integral_value = trapz(x, y);```三、线性方程组求解线性方程组求解是数值分析中的重要问题,常用的求解方法有高斯消元法和LU 分解法。
在Matlab中,可以使用\操作符来求解线性方程组。
例如,对于线性方程组Ax=b,可以使用如下代码进行求解:```matlabA = [1, 2; 3, 4];b = [5; 6];x = A\b;```四、实验结果与分析在本次实验中,我们分别使用Matlab进行了数值微分、数值积分和线性方程组求解的计算。
通过实验结果可以发现,Matlab提供了丰富的数值计算函数和工具,能够方便地进行数值分析的计算和求解。
数值微分的计算结果与解析解相比较,可以发现数值微分的误差随着步长的减小而减小,但是当步长过小时,数值微分的误差会受到舍入误差的影响。
matlab数值分析03-方程组求解
![matlab数值分析03-方程组求解](https://img.taocdn.com/s3/m/f8613ac5aa00b52acfc7cab4.png)
D = L = U = DL=
g = DL\b; tol = 1e-6; x = zeros(n,1); err = 1; while(err>tol) xt = DL\(U*x) + g; err = norm(xt-x); %keyboard; x = xt; end x
SOR迭代 SOR迭代
n = 10; A = rand(n)+ n*eye(n); b = A*ones(n,1); % For example, to solve Ax=b D = diag(diag(A)); L = -tril(A) + D; U = -triu(A) + D; omega = 0.6; DL= D-omega*L; DU= (1-omega)*D+omega*U; g = omega* (DL\b); tol = 1e-6; x = zeros(n,1); err = 1; while(err>tol) xt = DL\(DU*x) + g; err = norm(xt-x); %keyboard; x = xt; end
采用如下方式可求得精确的解: x1=inv(sym(A))*B x1 = [ -9/5, 12/5] [ 28/15, -19/15] [ 58/15, -49/15] [ -32/15, 41/15] norm(double(A*x1-B)) ans = 0
(2)若方程组有无穷多解 ) 齐次方程) Matlab语言中可以由null()直接求出, 语言中可以由null()直接求出 (齐次方程) 在Matlab语言中可以由null()直接求出, 其调用格式为: 其调用格式为: Z=null( 求解A的化0 Z=null(A) 求解A的化0矩阵 Z=null( Z=null(A,‘r’) 求解A的化0矩阵的规范形式 ) 求解A的化0 null()函数也可用于符号变量描述方程的解析解问题,其中Z ()函数也可用于符号变量描述方程的解析解问题 注:null()函数也可用于符号变量描述方程的解析解问题,其中Z 的列数为n 而各列构成的向量又称为矩阵A的基础解系。 的列数为n-r,而各列构成的向量又称为矩阵A的基础解系。 求解非齐次方程组也是比较简单的, 求解非齐次方程组也是比较简单的,只要能求出该方程的任 意一个特解x0 则原非齐次方程组的解为x=x1+x0 其实, x=x1+x0。 意一个特解x0 则原非齐次方程组的解为x=x1+x0。其实,在 Matlab中求解该方程的一个特解并非难事 中求解该方程的一个特解并非难事, Matlab中求解该方程的一个特解并非难事,用 x0=pinv( x0=pinv(A)*B 即可求出
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
佛山科学技术学院
实 验 报 告
课程名称 数值分析
实验项目 数值积分
专业班级 机械工程 姓 名 余红杰 学 号 2111505010 指导教师 陈剑 成 绩 日 期 月 日
一、实验目的
1、 掌握程序的录入和matlab 的使用和操作;
2、 了解影响线性方程组解的精度的因素——方法与问题的性态。
3、 学会Matlab 提供的“\”的求解线性方程组。
二、实验要求
1、按照题目要求完成实验内容;
2、写出相应的Matlab 程序;
3、给出实验结果(可以用表格展示实验结果);
4、分析和讨论实验结果并提出可能的优化实验。
5、写出实验报告。
三、实验步骤
1、用LU 分解及列主元高斯消去法解线性方程组
a)⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=⎪⎪⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎪⎪⎭⎫ ⎝⎛----15900001.582012151526099999.2310
7104321x x x x , 输出b Ax =中系数LU A =分解的矩阵L 和U ,解向量x 和)det(A ;用列主元法的行交换次序解向量x 和求)det(A ;比较两种方法所得结果。
2、用列主高斯消元法解线性方程组b Ax =。
(1)、⎪⎪⎪⎭
⎫ ⎝⎛=⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎭⎫ ⎝⎛--11134.981.4987.023.116.427
.199.103.601.3321x x x
(2)、⎪⎪⎪⎭
⎫ ⎝⎛=⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎭⎫ ⎝⎛--11134.981.4990.023.116.427
.199.103.600.3321x x x 分别输出)det(,,A b A ,解向量x ,(1)中A 的条件数。
分析比较(1)、(2)的计算结果
3、线性方程组b Ax =的A 和b 分别为
⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=1095791068565778710A ,⎪⎪⎪⎪⎪⎭
⎫ ⎝⎛=31332332b 则解T x ),1,1,1,1(=. 用MATLAB 内部函数求)det(A 和A 的所有特征值和2)(A cond . 若令
⎪⎪⎪⎪⎪⎭
⎫ ⎝⎛=+98.99599.6989.998.585604.508.72.71.8710A A δ, 求解b x x A A =++))((δδ,输出向量x δ和2x δ,从理论结果和实际计算两方面分析线性
方程组b Ax =解的相对误差22/x x
δ以及A 的相对误差
/A A δ的关系。
四、实验结果
1:
%run311.m
clc,clear;
A = [10 -7 0 1;-3 2.099999 6 2;5 -1 5 -1;2 1 0 2];
b = [8;5.90001;5;1];
%L U 分解
format short
%小数点后四位,不然会受到后面的影响
[L U] = lu(A)
%解方程组,输出A ,det(A)
y = L\b;
format long
%小数点后15位显示
x = U\y
%V ol_xiao.m列主元消去,未用增广矩阵而把系数矩阵A和矩阵B分开对应来的function x=V ol_xiao(A,b)
%列主元消去
%x为解
%A为系数矩阵
n =length(A);
%A默认为方阵,求其大小,且默认和b长度一样
C =zeros(1,n);
%用于系数矩阵交换存放数据
b0 =0;
%用于b矩阵交换存放数据
x =zeros(n,1);
%建一个矩阵存放解
%下面对矩阵进行n-1次列主元交换得到上三角矩阵
for i=1:(n-1)
for j=(i+1):n
if abs(A(i,i))<abs(A(j,i))
C =A(i,1:n);
A(i,1:n)=A(j,1:n);
A(j,1:n)=C;
b0 = b(i);
b(i) =b(j);
b(j) =b0;
else continue;
end
end;
for j=(i+1):n
b(j) = b(j) + b(i)*(-A(j,i))/A(i,i);
A(j,1:n)=A(j,1:n)+A(i,1:n).*(-A(j,i)/A(i,i));
end
end
%以下为判断是否满秩然后从下往上求解:
if A(n,n)~=0
x(n) = b(n)/A(n,n);
for i=(n-1):-1:1
sum = 0;
for j=(i+1):n
sum = sum + A(i,j)*x(j);
end
x(i) = (b(i)-sum)/A(i,i);
end
else x = ['err'];
end
%run312.m
clc,clear;
A = [10 -7 0 1;-3 2.099999 6 2;5 -1 5 -1;2 1 0 2];
b = [8;5.90001;5;1];
format long
x=V ol_xiao(A,b)
(其实可以看出,两种结果是一致的)
2.
%run321.m
clc;clear;
A = [3.01 6.03 1.99;1.27 4.16 -1.23;0.987 -4.81 9.34];
b = [1;1;1];
format bank
%格式转换为普通的,不然会用科学计数法表示
x = Vol_xiao(A,b)
T1 = cond(A,1)
T2 = cond(A,2)
T3 = cond(A,inf)
%可以看出条件数很大
%run322.m
clc;clear;
A = [3.01 6.03 1.99;1.27 4.16 -1.23;0.990 -4.81 9.34]
b = [1;1;1]
det1 = det(A)
format bank
x = Vol_xiao(A,b)
%从条件数就可以猜出其实微小变化就会对结果造成很大差异
%结果也验证了条件数大的方程的不稳定性很差。
3.
%run331.m
clc;clear;
A = [10 7 8 7;7 5 6 5;8 6 10 9;7 5 9 10];
b = [32;23;33;31];
format short
%显示对应行列式值和特征值
det1 = det(A)
eig1 = eig(A)
format bank
%显示条件数
tjs2 = cond(A,2)
%run332.m
clc;clear;
A = [10 7 8 7;7 5 6 5;8 6 10 9;7 5 9 10];
b = [32;23;33;31];
X = [1;1;1;1];
AA= [10 7 8.1 7.2;7.08 5.04 6 5;8 5.98 9.89 9;6.99 5 9 9.98];
oA= AA - A;
XX= AA \ b;
oX= XX – X
%输出oX 的2范数
normx2 = norm(oX)
%输出X 的相对误差
xrerro = normx2/norm(X)
%输出A 的相对误差
Arerro = norm(oA)/norm(A)
五、讨论分析
1.对于第一题,有一个顺序主子式是极小的,接近等于0,所以无法直接得到标准的上下三角,所以调用程序内置的lu 函数进行分解,所求到的值与列主元方法的值一致,查找后发现matlab 的lu 分解就是根据列主元方法求出的;
2.
条件数对一个方程的稳定性有很大影响,题中的系数矩阵条件数超过五万,微
小的变化就引起了结果的很大变化,如果是数据分析的话相当于蝴蝶效应般严重;
3.第三题也是对于一个系数矩阵的范数,很小的变化引起结果的大变化,该系数矩阵的相关数也是三千多,也是比较大,不稳定的。
六、改进实验建议
相关题目太少,而且很少有扩展性题目,对于各种指令的使用并不充分,其实在这个求解方程组的过程中,还是有比较多的相关方法可以进行MATLAB实现的。