MATLAB中解三对角线性方程组追赶法程序

合集下载

(完整word版)线性方程组的直接解法及matlab的实现

(完整word版)线性方程组的直接解法及matlab的实现

本科毕业论文(2010 届)题目线性方程组的直接解法及matlab的实现学院数学与信息工程学院专业数学与应用数学班级2006级数学1 班学号0604010127学生姓名胡婷婷指导教师王洁完成日期2010年5月摘要随着科技技术的发展及人类对自然界的不断探索模拟。

在自然科学和工程问题中的很多问题的解决常常归结为线性代数问题!本文的主要内容是对线性方程组求解方法的探讨,主要介绍了四种求解线性方程组的方法,第一种是教科书上常见的消元法,我们称之为基本法。

第二种方法是标准上三角形求解法,即将增广矩阵经过初等变换后化成标准上三角形,然后求解.它改进了一般教科书上的常见方法,与常见方法比较有如下优点:1)规范了自由未知量的选取;2)只用矩阵运算;3)减少了计算量.第三种方法是对特定的方程组(系数矩阵A为n阶对称正定矩阵,且A的顺序主子式均不为零。

)的求解方法进行描述,并且为这种线性方程的求解提供了固定的公式化的方法。

第四种方法是对现在实际问题中常常会遇到的系数矩阵为三对角矩阵的方程组的求解方法。

同时给出这几种方法的数值解法(matlab程序),由于运用电脑软件求解,所以必须考虑计算方法的时间、空间上的效率以及算法的数值稳定性问题,所以针对不同类型的线性方程组有不同的解法.但是,基本的方法可以归结为两大类,即直接法和迭代法.关键词高斯消去法;三角分解法;乔莱斯基分解法;追赶法AbstractSystems of linear equations are associated with many problems in engineering and scinence ,as well as with applications of mathematics to the social sciences and the quantitative study of business and economic problems.The main content of this article is the method for solving linear equations,we introduce four methods for solving linear equations in this paper。

matlab追赶法解常微分方程

matlab追赶法解常微分方程

研究领域:数学、计算机科学文章标题:深入探讨matlab追赶法解常微分方程在数学和计算机科学领域中,常微分方程是一个重要且广泛应用的课题。

而matlab追赶法作为常微分方程的求解方法,在实际应用中具有重要意义。

本文将以深度和广度兼具的方式,对matlab追赶法解常微分方程这一主题展开全面评估,并撰写一篇有价值的文章,同时结合个人观点和理解,为读者提供深刻的思考。

一、matlab追赶法解常微分方程简介1.1 matlab追赶法基本原理matlab追赶法,又称托马斯算法,是一种用于求解三对角线性方程组的方法。

在常微分方程的数值解法中,常常会遇到需要求解三对角线性方程组的情况,而matlab追赶法正是针对这一问题而提出的高效算法。

1.2 追赶法在常微分方程求解中的应用常微分方程在实际问题中有着广泛的应用,而求解常微分方程的过程中往往需要用到追赶法。

追赶法不仅可以提高计算效率,还可以有效地解决数值稳定性和精度的问题,因此在工程和科学计算中得到了广泛的应用。

二、深入探讨matlab追赶法解常微分方程2.1 算法实现及优化matlab追赶法的实现涉及到矩阵运算、追赶过程和追赶系数的求解等关键步骤。

如何针对不同类型的方程组进行算法优化,是一个需要深入探讨的问题。

通过优化算法,可以提高追赶法的计算效率和数值稳定性,使其在常微分方程求解中发挥更大的作用。

2.2 算法的数值分析通过数值分析,可以更加深入地了解matlab追赶法在解常微分方程过程中的数值特性。

包括收敛性、稳定性、误差分析等方面,这些都是影响算法性能和应用效果的重要因素,需要进行深入的研究和分析。

三、对matlab追赶法解常微分方程的个人观点和理解3.1 算法的优势与局限性matlab追赶法作为一种高效的求解算法,具有较好的稳定性和精度,特别适合于大规模的常微分方程求解。

但在某些特定问题上,追赶法的适用性和效率仍然存在局限性,需要进行合理的选择和应用。

数值分析中求解线性方程组的MATLAB程序(6种)

数值分析中求解线性方程组的MATLAB程序(6种)

数值分析中求解线性方程组的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作业

#数值分析Matlab作业

数值分析编程作业2012年12月第二章14.考虑梯形电阻电路的设计,电路如下:电路中的各个电流{i1,i2,…,i8}须满足下列线性方程组:121232343454565676787822/252025202520252025202520250i i V R i i i i i i i i i i i i i i i i i i i i -=-+-=-+-=-+-=-+-=-+-=-+-=-+=这是一个三对角方程组。

设V=220V ,R=27Ω,运用追赶法,求各段电路的电流量。

Matlab 程序如下:function chase () %追赶法求梯形电路中各段的电流量 a=input('请输入下主对角线向量a='); b=input('请输入主对角线向量b='); c=input('请输入上主对角线向量c='); d=input('请输入右端向量d='); n=input('请输入系数矩阵维数n='); u(1)=b(1); for i=2:nl(i)=a(i)/u(i-1); u(i)=b(i)-c(i-1)*l(i); endy(1)=d(1); for i=2:ny(i)=d(i)-l(i)*y(i-1); endx(n)=y(n)/u(n); i=n-1; while i>0x(i)=(y(i)-c(i)*x(i+1))/u(i); i=i-1;end x输入如下: >> chase请输入下主对角线向量a=[0,-2,-2,-2,-2,-2,-2,-2]; 请输入主对角线向量b=[2,5,5,5,5,5,5,5];请输入上主对角线向量c=[-2,-2,-2,-2,-2,-2,-2,0]; 请输入方程组右端向量d=[220/27,0,0,0,0,0,0,0]; 请输入系数矩阵阶数n=8 运行结果如下:x = 8.1478 4.0737 2.0365 1.0175 0.5073 0.2506 0.1194 0.0477第三章14.试分别用(1)Jacobi 迭代法;(2)Gauss-Seidel 迭代法解线性方程组1234510123412191232721735143231211743511512x x x x x ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥---⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥=--⎢⎥⎢⎥⎢⎥--⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥---⎣⎦⎣⎦⎣⎦ 迭代初始向量(0)(0,0,0,0,0)T x =。

matlab线性方程组的矩阵解法

matlab线性方程组的矩阵解法
2 解 Ux = −1 7 2 0 0 1 -1 x1 1 x1 13 1 0 x 2 2 x2 2 ,得 2 = 13 x = − 13 ,得 x 0 −1 3 7 3 7 x x4 −1 0 11 4 − 11 7 7 0
function x=lupdsv(A,b) n=length(b); [LU,p]=lupd(A); y(1)=b(p(1)); for i=2:n y(i)=b(p(i))-LU(i,1:i-1)*y(1:i-1)'; end x(n)=y(n)/LU(n,n); for i=(n-1):-1:1 x(i)=(y(i)-LU(i,i+1:n)*x(i+1:n)')/LU(i,i); end
lupdsv.m %功能:调用列主元三角分解函数 [LU,p]=lupd(A) % 求解线性方程组Ax=b。 。 求解线性方程组
%解法:PA=LU, Ax=b←→PAx=Pb 解法: 解法 % % LUx=Pb, Ly=f=Pb, y=Ux f(i)=b(p(i))
%输入:方阵A,右端项 (行或列向量均可) 输入:方阵 ,右端项b(行或列向量均可) 输入 %输出:解x(行向量) 输出: 输出 (行向量)
Ax = d 用矩阵表示 应用追赶法求解三对角线性方程组。追赶法仍然 追赶法求解三对角线性方程组 应用追赶法求解三对角线性方程组。追赶法仍然 保持LU分解特性,它是一种特殊的LU分解。 LU分解特性 LU分解 保持LU分解特性,它是一种特殊的LU分解。充分利用 了系数矩阵的特点,而且使之分解更简单, 了系数矩阵的特点,而且使之分解更简单,得到对三对 角线性方程组的快速解法。 角线性方程组的快速解法。

matlab求解块三对角矩阵

matlab求解块三对角矩阵

matlab求解块三对角矩阵求解块三对角矩阵是数值线性代数领域中的一个经典问题。

在本文中,我们将使用MATLAB来解决这个问题,并详细介绍解决方法和步骤。

块三对角矩阵是一种特殊的矩阵结构,其中除了主对角线上的元素外,还存在两个对角线上的块状子矩阵。

这种矩阵在科学和工程计算中经常出现,如有限元方法、微分方程数值解等。

我们需要定义一个具体的块三对角矩阵。

假设我们的矩阵是n×n的矩阵,其中每个块状子矩阵是m×m的矩阵。

为了简化问题,我们假设矩阵的每个子矩阵都是对角矩阵,即只有主对角线上有非零元素。

接下来,我们需要将块三对角矩阵表示为MATLAB中的矩阵形式。

可以使用MATLAB中的稀疏矩阵来表示这种特殊的矩阵结构。

稀疏矩阵是一种只存储非零元素的矩阵表示方法,可以大大减少存储空间和计算时间。

在MATLAB中,可以使用spdiags函数来创建块三对角矩阵。

该函数需要两个参数:一个是主对角线上的元素,另一个是两个对角线上的块状子矩阵。

我们可以通过定义这两个参数来创建我们需要的块三对角矩阵。

接下来,我们需要解决块三对角矩阵的求解问题。

求解块三对角矩阵的方法有很多种,如追赶法、LU分解法、雅可比迭代法等。

在本文中,我们将使用追赶法来求解块三对角矩阵。

追赶法是一种特定的高斯消元法,适用于求解具有特殊结构的线性方程组。

它通过对矩阵进行分解,将原问题转化为两个简单的三对角线性方程组的求解问题。

在MATLAB中,可以使用函数cholfact来对块三对角矩阵进行分解。

该函数返回一个分解对象,可以用于求解线性方程组。

然后,我们可以使用该对象的solve方法来求解线性方程组。

在实际求解过程中,我们需要定义一个右侧向量b,并将其传递给solve方法。

solve方法将返回一个解向量x,满足Ax=b,其中A是我们的块三对角矩阵。

我们可以将解向量x打印出来,并进行验证。

可以通过计算残差向量r=Ax-b来验证求解的准确性。

matlab追赶法解101阶三对角方程组

matlab追赶法解101阶三对角方程组

在探讨MATLAB追赶法解101阶三对角方程组之前,我们首先需要了解什么是追赶法和什么是三对角方程组。

追赶法又称托马斯算法,是一种用于求解带状矩阵(即只有主对角线和两条相邻的对角线上有非零元素的矩阵)的线性方程组的方法。

而三对角矩阵就是只有主对角线和两条相邻的对角线上有非零元素的矩阵。

在实际应用中,求解带状矩阵的线性方程组是非常常见的,特别是在数值计算和科学工程领域。

现在,让我们深入探讨MATLAB追赶法解101阶三对角方程组的方法和具体步骤。

一、MATLAB追赶法解101阶三对角方程组1. 概念介绍101阶三对角方程组是一个非常大的线性方程组,通常使用传统的高斯消元法来求解会耗费大量的时间和计算资源。

而MATLAB追赶法通过利用三对角矩阵的特殊性质,可以有效地简化计算过程,并且节省大量的内存和计算资源。

2. 追赶法步骤(1)将原方程组化为追赶法所需的形式;(2)利用追赶法求解三对角线性方程组。

二、追赶法求解101阶三对角方程组的实现过程1. 将原方程组化为追赶法所需的形式对于101阶三对角方程组,我们首先需要将其化为追赶法所需的形式。

这个过程涉及到选取合适的追赶元和追赶子以及对原方程组的变形,将其化为追赶法能够直接处理的形式。

2. 利用追赶法求解线性方程组一旦将原方程组化为追赶法所需的形式,我们就可以利用追赶法对其进行求解。

追赶法的核心是通过追赶子的迭代计算,逐步求得线性方程组的解。

在MATLAB中,可以使用内置的追赶法求解函数,也可以编写自定义的追赶法算法来实现对101阶三对角方程组的求解。

三、个人观点和理解在实际工程和科学计算中,追赶法是一种非常有效的求解带状矩阵线性方程组的方法。

对于大规模的三对角方程组,特别是高阶的情况,传统的直接求解方法往往会遇到内存和计算资源的限制,而追赶法能够通过精巧的迭代计算,在保证解的精度的显著提高计算效率。

在MATLAB中,通过调用内置的追赶法函数,可以快速地求解大规模的三对角方程组,极大地方便了工程实践中的数值计算工作。

用matlab解线性方程组

用matlab解线性方程组

用matlab解线性方程组2008-04-12 17:00一。

高斯消去法1.顺序高斯消去法直接编写命令文件a=[]d=[]'[n,n]=size(a);c=n+1a(:,c)=d;for k=1:n-1a(k+1:n, k:c)=a(k+1:n, k:c)-(a(k+1:n,k)/ a(k,k))*a(k, k:c); %消去endx=[0 0 0 0]' %回带x(n)=a(n,c)/a(n,n);for g=n-1:-1:1x(g)=(a(g,c)-a(g,g+1:n)*x(g+1:n))/a(g,g)end2.列主高斯消去法*由于“[r,m]=max(abs(a(k:n,k)))”返回的行是“k:n,k”内的第几行,所以要通过修正来把m 改成真正的行的值。

该程序只是演示程序,真正机器计算不需要算主元素所在列以下各行应为零的值。

直接编写命令文件a=[]d=[] '[n,n]=size(a);c=n+1a(:,c)=d; %(增广)for k=1:n-1[r,m]=max(abs(a(k:n,k))); %选主m=m+k-1; %(修正操作行的值)if(a(m,k)~=0)if(m~=k)a([k m],:)=a([m k],:); %换行enda(k+1:n, k:c)=a(k+1:n, k:c)-(a(k+1:n,k)/ a(k,k))*a(k, k:c); %消去endendx=[0 0 0 0]' %回带x(n)=a(n,c)/a(n,n);for g=n-1:-1:1x(g)=(a(g,c)-a(g,g+1:n)*x(g+1:n))/a(g,g)end3.分别用顺序高斯消去法和列主高斯消去法解方程组a*x=d,并比较结果a=[0 1 2 3;9 11 23 34;62.5 23.4 15.5 17.2;192.01 124 25.1 59.3] d=[1;1;1;1]顺序高斯消去法:提示“Warning: Divide by zero.” x =NaN NaN NaN NaN 列主高斯消去法:x =-1.2460 2.8796 5.5206 -4.3069由此可见列主高斯消去法可以解决顺序高斯消去法所不能解决的问题。

追赶法求解三对角线性方程组

追赶法求解三对角线性方程组

追赶法求解三对角线性方程组一 实验目的利用编程方法实现追赶法求解三对角线性方程组。

二 实验内容1、 学习和理解追赶法求解三对角线性方程组的原理及方法;2、 利用MA TLAB 编程实现追赶法;3、 举例进行求解,并对结果进行分。

三 实验原理设n 元线性方程组Ax=d 的系数矩阵A 为非奇异的三对角矩阵11222=(1)(n 1)()()a c b a c A a n c b n a n ⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥--⎢⎥⎢⎥⎣⎦………… 这种方程组称为三对角线性方程组。

显然,A 是上下半宽带都是1的带状矩阵。

设A 的前n-1个顺序主子式都不为零,根据定理2.5的推论,A 有唯一的Crout 分解,并且是保留带宽的。

其中L 是下三角矩阵,U 是单位上三角矩阵。

利用矩阵相乘法,可以1112212(1)1u(n 1)()()1l u m l u A LU l n m n l n ⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥==⨯⎢⎥⎢⎥--⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦……………得到:由上列各式可以得到L 和U 。

引入中间量y ,令yUx =,则有:已知L 和d ,可求得y 。

则可得到y 的求解表达式:11/12,3,,()(1)*y()=()[()(1)]/y d l i nm i y i li i di y i di m i y i li==-+=--…1111111/1(2)(1)(1)u (1)(11)/(1)(1)(1)l a l u c u c l mi bi i n a i m i i l i i n ci li ui ui ci li l i a i b i ui=*===≤≤+=+++≤≤-=∙=+=+-+Ax LUx Ly d Ly d ====1112222(1)(n 1)(n 1)()()(n)(n)l y d m l y d l n y d m n l n y d ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⨯=⎢⎥⎢⎥⎢⎥---⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦……………由y Ux =得:111112221u(n 1)(n 1)(n 1)1(n)(n)u x y u x y x y x y ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⨯=⎢⎥⎢⎥⎢⎥---⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦………… 可得到X 的求解表达式:()()1,2,,1()()u()(1)x n y n i n n x i y i i x i ==--=-+… 从而得到Ax=d 的解x 。

MATLAB-追赶法求解三对角方程组的算法原理例题与程序

MATLAB-追赶法求解三对角方程组的算法原理例题与程序

3)三对角形线性方程组123456789104100000000141000000001410000000014100000000141000000001410000000014100000000141000000001410000000014x x x x x x x x x x -⎡⎤⎡⎤⎢⎥⎢⎥--⎢⎥⎢⎥⎢⎥⎢⎥--⎢⎥⎢⎥--⎢⎥⎢⎥⎢⎥⎢⎥--⎢⎥⎢⎥--⎢⎥⎢⎥⎢⎥⎢⎥--⎢⎥⎢⎥--⎢⎥⎢⎥⎢⎢⎥--⎢⎢⎥⎢⎢⎥-⎣⎦⎣⎦7513261214455⎡⎤⎢⎥⎢⎥⎢⎥-⎢⎥⎢⎥⎢⎥=⎢⎥-⎢⎥⎢⎥⎢⎥-⎢⎥⎥⎢⎥⎥⎢⎥⎥⎢⎥-⎣⎦*(2,1,3,0,1,2,3,0,1,1)T x =--- 二、数学原理设系数矩阵为三对角矩阵112223311100000000000000000n n n n n b c a b c a b A a b c a b ---⎛⎫ ⎪⎪ ⎪=⎪ ⎪ ⎪⎪⎪⎝⎭L L L M M MM M M L L则方程组Ax=f 称为三对角方程组。

设矩阵A 非奇异,A 有Crout 分解A=LU ,其中L 为下三角矩阵,U 为单位上三角矩阵,记11222331100001000000010000000100,00000000000001n n nn b L U γαβγββγβ--⎛⎫⎛⎫⎪ ⎪⎪ ⎪ ⎪ ⎪∂==⎪⎪ ⎪ ⎪ ⎪ ⎪⎪ ⎪⎪ ⎪∂⎝⎭⎝⎭L L L LL L M M M MM M M MM M L L LLL可先依次求出L ,U 中的元素后,令Ux=y ,先求解下三角方程组Ly=f 得出y ,再求解上三角方程组Ux=y 。

事实上,求解三对角方程组的2追赶法将矩阵三角分解的计算与求解两个三角方程组的计算放在一起,使算法更为紧凑。

其计算公式为:1111,1111,111,2,3,,,1,2,,1ii i i i i i i ii i i i i n ni i i i c f b y i n c a b a f y y x y i n n x y x βγββαβγγβαβγ--+⎧===⎪⎪=⎪⎪⎪==-=⎪⎪⎨-⎪=⎪⎪=⎪⎪=--⎪=-⎪⎩L L 对对(*)三、程序设计function x=chase(a,b,c,f)%求解线性方程组Ax=f,其中A 是三对角阵 %a 是矩阵A 的下对角线元素a(1)=0 %b 是矩阵A 的对角线元素%c 是矩阵A 的上对角线元素c(n)=0 %f 是方程组的右端向量 n=length(f);x=zeros(1,n);y=zeros(1,n); d=zeros(1,n);u= zeros(1,n); %预处理 d(1)=b(1); for i=1:n-1 u(i)=c(i)/d(i);d(i+1)=b(i+1)-a(i+1)*u(i); end%追的过程y(1)=f(1)/d(1); for i=2:ny(i)=(f(i)-a(i)*y(i-1))/d(i); end%赶的过程 x(n)=y(n); for i=n-1:-1:1x(i)=y(i)-u(i)*x(i+1); end>> a=[0,-1,-1,-1,-1,-1,-1,-1,-1,-1];>> b=[4,4,4,4,4,4,4,4,4,4];>> c=[-1,-1,-1,-1,-1,-1,-1,-1,-1,0];>> f=[7,5,-13,2,6,-12,14,-4,5,-5];>> x=chase(a,b,c,f)x =2.00001.0000-3.00000.00001.0000-2.00003.0000-0.00001.0000-1.0000四、结果分析和讨论追赶法求解的结果为x=(2,1,-3,0,1,-2,3,0,1,-1)T。

matlab线性方程组的矩阵解法

matlab线性方程组的矩阵解法
用矩阵分解法求解线性方程组
一、利用三角分解A = LU 求解Ax = b
二、 用列主元的三角分解PA = LU 求解Ax = b
三、 用全主元的三角分解PAQT = LU 求解Ax = b
四、利用Cholesky分解A = LLT 求解Ax = b
五、 三对角方程组的解法
一、利用三角分解A = LU 求解Ax = b 设已有 A = LU 代入原方程 Ax = b 得 LY = b LUx = b ⇔ Ux = Y
定理 如果带宽为 w=p+q-1 的n阶带状矩阵 有LU 阶带状矩阵A有 阶带状矩阵 分解:A=LU,则L是带宽为 的下三角矩阵,U是带宽 是带宽为p的下三角矩阵 分解 , 是带宽为 的下三角矩阵, 是带宽 的上三角矩阵。 为q的上三角矩阵。 的上三角矩阵
a11 M a p1 0 L a22 O O an , n − p + 1 a1q O an− q +1,n O M L ann 0 L u22 u1q O un − q +1,n O O M unn 0
j = k +1
∑u
n
kj
x j ) ukk y1 y 2 M yn
( k = n − 1, n − 2, L ,1) u11 u12 L u1 n x1 u22 L u2 n x 2 = O M M unn x n
Ax = d 用矩阵表示 应用追赶法求解三对角线性方程组。追赶法仍然 追赶法求解三对角线性方程组 应用追赶法求解三对角线性方程组。追赶法仍然 保持LU分解特性,它是一种特殊的LU分解。 LU分解特性 LU分解 保持LU分解特性,它是一种特殊的LU分解。充分利用 了系数矩阵的特点,而且使之分解更简单, 了系数矩阵的特点,而且使之分解更简单,得到对三对 角线性方程组的快速解法。 角线性方程组的快速解法。

高斯消去、追赶法matlab

高斯消去、追赶法matlab

⾼斯消去、追赶法matlab 1. 分别⽤Gauss消去法、列主元Gauss消去法、三⾓分解⽅法求解⽅程组程序:(1)Guess消去法:function x=GaussXQByOrder(A,b)%Gauss消去法N = size(A);n = N(1);x = zeros(n,1);for i=1:(n-1)for j=(i+1):nif(A(i,i)==0)disp('对⾓元不能为0');return;endm = A(j,i)/A(i,i);A(j,i:n)=A(j,i:n)-m*A(i,i:n);b(j)=b(j)-m*b(i);endendx(n)=b(n)/A(n,n);for i=n-1:-1:1x(i)=(b(i)-sum(A(i,i+1:n)*x(i+1:n)))/A(i,i);end命令⾏输⼊:A=[1 -1 2 1;-1 3 0 -3 ;2 0 9 -6;1 -3 -6 19];b=[1 3 5 7];b=b';x=GaussXQByOrder(A,b)运算结果:x =-8.0000000000000000.3333333333333333.6666666666666672.000000000000000(2)列主元Gauss消去法程序:function x=GaussXQLineMain(A,b)%列主元Gauss消去法N = size(A);n = N(1);x = zeros(n,1);zz=zeros(1,n);for i=1:(n-1)[~,p]=max(abs(A(i:n,i)));zz=A(i,:);A(i,:)=A(p+i-1,:);A(p+i-1,:)=zz;temp=b(i);b(i)=b(i+p-1);b(i+p-1)=temp;for j=(i+1):nm = A(j,i)/A(i,i);A(j,i:n)=A(j,i:n)-m*A(i,i:n);b(j)=b(j)-m*b(i);endendx(n)=b(n)/A(n,n);for i=n-1:-1:1x(i)=(b(i)-sum(A(i,i+1:n)*x(i+1:n)))/A(i,i);end命令⾏:A=[1 -1 2 1;-1 3 0 -3 ;2 0 9 -6;1 -3 -6 19];b=[1 3 5 7];b=b';x=GaussXQLineMain(A,b)运⾏结果:x =-8.0000000000000050.3333333333333323.6666666666666682.000000000000000(3)三⾓分解⽅法程序:function x = LU(A,b)%三⾓分解N = size(A);n = N(1);L = eye(n,n);U = zeros(n,n);x = zeros(n,1);y = zeros(n,1);U(1,1:n) = A(1,1:n);L(1:n,1) = A(1:n,1)/U(1,1);for k=2:nfor i=k:nU(k,i) = A(k,i)-L(k,1:(k-1))*U(1:(k-1),i);endfor j=(k+1):nL(j,k) = (A(j,k)-L(j,1:(k-1))*U(1:(k-1),k))/U(k,k); endendy(1)=b(1)/L(1,1);for i=2:ny(i)=b(i)-sum(L(i,1:i-1)*y(1:i-1));endx(n)=y(n)/U(n,n);for i=n-1:-1:1x(i)=(y(i)-sum(U(i,i+1:n)*x(i+1:n)))/U(i,i);end命令⾏:A=[1 -1 2 1;-1 3 0 -3 ;2 0 9 -6;1 -3 -6 19];b=[1 3 5 7];b=b';x=LU(A,b)运⾏结果:x =-8.0000000000000000.3333333333333333.6666666666666672.000000000000000程序:function [times,wucha]=zhuiganfa(a,b,c,f)%追赶法:x为所求解,times为所有乘除运算次数(即时间),wucha为误差的2-范数。

matlab中调用追赶法求解系数矩阵的三次样条插值

matlab中调用追赶法求解系数矩阵的三次样条插值

matlab中调用追赶法求解系数矩阵的三次样条插值已知数表x1245y1342边界条件S''(x0)=0,S''(x3)=0插值点x=3---------------------------------------------------------分--割--线---------------------------------------------------------使用三转角方程边界条件为第二类边界条件---------------------------------------------------------分--割--线---------------------------------------------------------先写好追赶法的程序把followup.m存入工作路径function x=followup(a,b,c,d)n=length(d);a(1)=0;%“追”的过程L(1)=b(1);y(1)=d(1)/L(1);u(1)=c(1)/L(1);for i=2:(n-1)L(i)=b(i)-a(i)*u(i-1);y(i)=(d(i)-y(i-1)*a(i))/L(i);u(i)=c(i)/L(i);endL(n)=b(n)-a(n)*u(n-1);y(n)=(d(n)-y(n-1)*a(n))/L(n);%“赶”的过程x(n)=y(n);for i=(n-1):-1:1x(i)=y(i)-u(i)*x(i+1);end---------------------分--割--线用matlab编程如下:function[s,y0]=spline3 (x,y,x0)%x,y为数表x0为插值点s表示插值函数y0为x0对应的插值函数值syms tn=length(x);%得出nfor i=1:n-1;h(i)=x(i+1)-x(i);endfor i=2:n-1;lamda(i)=h(i)/(h(i-1)+h(i));miu(i)=1-lamda(i);g(i)=3*(lamda(i)*((y(i)-y(i-1))/h(i-1))+miu(i)*((y(i+1)-y(i))/h( i)));endg(1)=3*((y(2)-y(1))/h(1));g(n)=3*((y(n)-y(n-1))/h(n-1));%前边求出lamda miu和g从而可以确定系数矩阵miu(1)=1;miu(4)=0;lamda(n)=1;%根据第二边界条件补充两个lamda和miu的值for i=1:nbeta(i)=2;endm=followup(lamda,beta,miu,g)%解出m的值从而可确定st st为各段的插值多项式for i=1:n-1st(i)=(t-x(i+1))^2*(h(i)+2*(t-x(i)))*y(i)/(h(i)^3)...+(t-x(i))^2*(h(i)+2*(x(i+1)-t))*y(i+1)/(h(i)^3)...+(t-x(i))^2*(t-x(i+1))*m(i+1)/(h(i)^2)...+(t-x(i+1))^2*(t-x(i))*m(i)/(h(i)^2);end%得到插值的结果各段的t的表达式%接下来要将插值点x0代入首先确定x0所在的插值区间for i=1:n-1if (x(i)<x0)&&(x(i+1)>x0)in=i;endends=st(in);s=expand(s);s=collect(s,'t');y0=subs(s,'t',x0)%s是插值多项式y0是插值点的函数值---------------------------------------------------------分--割--线---------------------------------------------------------在matlab中输入x=[1 2 4 5];y=[1 3 4 2];spline3(x,y,3)---------------------------------------------------------分--割--线---------------------------------------------------------得到y0 =4.2500ans =-1/8*t^3+3/8*t^2+7/4*t-1。

MATLAB计算方法3解线性方程组计算解法名师公开课获奖课件百校联赛一等奖课件

MATLAB计算方法3解线性方程组计算解法名师公开课获奖课件百校联赛一等奖课件

li1 ai1
u11
(i 2,3,, n)
k 1
ukj akj lkmumj akj
m 1
(
j
k,
k
1,,
n)
lik
aik
k 1
limumk
m 1
(i
k
1,,
n)
ukk aik
(k 2,3,, n)
例3.1
2 1 2 6 2 1 2 6
4 5 4 18 2 3 0 6
a11 a12 a1n l11
a21
a22
a2n
l21
l22
l11 l21 l n1
l22
l
n2
an1
an2
ann
l n1
l n2
l
nn
l
nn
其中aij a ji
由矩阵乘法
(1)
1)
l2 11
a11
l11
a11
(取正)
2) L第1行 LT第j列 (j 2,,n)
…….
(k)
1求u的第k行:用L的第k行 u的第j列
(j k,k 1,,n)
(lk1 , lk 2 ,, lkk,0,0) (u1 j , u2 j ,, u jj,0,0)' akj
k 1
k 1
lkmumj 1 ukj akj ukj akj lkmumj
m 1
m 1
2 求L的第k列:用L的第i行 u的第k列
利用Gauss消元法得到同解旳三角方程为
1 c1
y1
2 c2
y2
n1
ቤተ መጻሕፍቲ ባይዱ
cn1

三对角方程组求解的追赶法

三对角方程组求解的追赶法

1
2 5
1
1
2 2
1 1
3
5
2
A2
1 2
0
1 2
1 2
0
1 2
例题
(3) k 3,z(3) 3,c3 2
1 1 3 4 1 3
A2 112
1 0
252
3 1
1 0
5 2 2
4 1 6
3 1
5
A3
1 0 2
例题
由z(2) 3,z(1) 3
1 1 1 0 1 0 0 1 3 1 1 0
2 2 1 0 0 1 0 2 5 2 0 1
例题
1 0 4 1 2 0 1 0 0 1 6 4
0
1
3 1
1
0
0
1
01
5
3
0 0 1 0 2 1 0 0 1 0 2 1
所以得
1 6 4
A1
1
5 3
0 2 1
解法二:Gauss Jordan原地求逆法
矩阵求逆
实现 :由Gauss Jordan消去法矩阵形式
Ln Ln1...L1A I
1
...
l(k)
1k
其中LK
l(k) kk
l(k) nk
1
故A1 Ln Ln1...L1
l (k )
ik
a(k) ik
a(k) kk
(i
k)
1 (i k)
a(k) kk
Gauss—Jordan原地求逆法
2
p3 1
y3
f3
pn 1 yn fn
解得
y1 yi
f1 fi

追赶法解三对角矩阵

追赶法解三对角矩阵

实验追赶法解三对角方程组一、实验目的学会用追赶法解三对角方程组,并应用该算法于实际问题.二、实验要求给定三对角方程组,应用追赶法解得方程组的解。

三、实验内容1、追赶法2、以课本数值试验2为实例3、如果有错,修改直至运行成功,查看运行结果;四、实验环境matlab五、实验步骤和方法1、程序设计2、带入实例3、撰写实验报告。

六、实验预习要求得到实例的解一、[源程序]function x = my_zgf2(A,d,flag)%MY_ZGF2 Summary of this function goes here[m,n]=size(A); %计算矩阵的大小if nargin==2; %输入变量等于2的时候,A中储存所有元素的值for i=1:na(i)=A(i+1,i);b(i)=A(i,i);c(i)=A(i,i+1);enda(1)=0; %补充不足的值b(n)=A(n,n);c(n)=0;elsec=[A(1,:) 0]; %flag==1时b=A(2,:);a=[0 A(3,:)];endu(1)=b(1);for i=2:n %第一次追赶,得到上、下三角矩阵l(i)=a(i)/u(i-1);u(i)=b(i)-c(i-1)*l(i);endy(1)=d(1); %解Ly=dfor i=2:ny(i)=d(i)-l(i)*y(i-1);endx(n)=y(n)/u(n); %解Ux=yfor i=n-1:-1:1x(i)=(y(i)-c(i)*x(i+1))/u(i);end二、带入实例A =-2.0000 -2.0000 -2.0000 -2.0000 -2.0000 -2.0000 -2.0000 05.0000 5.0000 5.0000 5.0000 5.0000 5.0000 5.0000 5.0000-2.0000 -2.0000 -2.0000 -2.0000 -2.0000 -2.0000 -2.0000 0d= 8.1400 0 0 0 0 0 0 0>> d=A(4,:);my_zgf2(A,d,1)ans =2.0350 1.0174 0.5086 0.2541 0.1267 0.0626 0.0298 0.0119 >>。

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