列主元高斯消去法和列主元三角分解法解线性方程
用列主元高斯消元法求线性代数方程组的解

课程设计任务书前 言回顾普通解方程组的方法,一般都是先逐个削去未知变量,最终得到只有一个未知变量的方程,解之,把得到的值回代到消去变量过程中得到的方程组,逐个求出未知变量。
这种解线性方程组的基本方法就是这里要介绍的高斯消去法。
数学上,高斯消元法(或译:高斯消去法),是线性代数中的一个算法,可用来为线性方程组求解,求出矩阵的秩,以及求出可逆方阵的逆矩阵。
当用于一个矩阵时,高斯消元法会产生出一个“行梯阵式”。
高斯消元法可以用在电脑中来解决数千条等式及未知数。
高斯消元法可以用来找出一个可逆矩阵的逆矩阵。
用关联矩阵表述网络拓扑结构,并根据厂站拓扑结构和网络拓扑结构等概念简化了电力系统的拓扑结构。
根据广义乘法和广义加法的运算规则,将改进的高斯消元算法应用于电力系统拓扑结构分析中,并引入稀疏、分块处理等技术提高了上述拓扑分析的效率。
采用上述高斯消元算法对山东电网220kV 以上的变电站进行拓扑结构分析,结果表明了运用该高斯消元法进行网络拓扑分析的正确性和有效性。
用列主元素法,选取每列的绝对值最大的元素作为消去对象并作为主元素。
然后换行使之变到主元位子上,在进行消元计算。
设)()(k k b X A ,确定第k 列主元所在位置k i ,在交换k i 行和k 行后,在进行消元,并用MATLAB 软件进行求解。
目录摘要....................................................................................... 错误!未定义书签。
第1章绪论 ......................................................................... 错误!未定义书签。
第2章高斯消元法的算法描述 (2)2.1高斯消元法的原理概述 (2)c231730658" 2.1.1高斯消元法的消元过程 (2)c231730658" 2.1.2高斯消元法的回带过程 (3)c231730658" 2.1.3高斯消元法的复杂度分析 (4)c231730658" 2.2列主高斯消元法原理简介 (5)c231730658" 2.2.1列主高斯消元法的消元过程 (6)c231730658" 2.2.2列主高斯消元法的回带过程 (6)c231730658" 2.2.3列主高斯消元法的算法描述 (6)c231730662"第3章高斯消元法的物理应用 (9)3.1c231730663"电网模型的描述 (9)c231730658" 3.2电网模型的问题分析 (9)c231730658"3.3求解计算 (11)c231730693"参考文献 (13)摘 要用列主元素高斯消去法法,选取每列的绝对值最大的元素作为消去对象并作为主元素。
多元变量的方程组求解

多元变量的方程组求解在许多实际问题中,常常需要求解由多个变量组成的方程组。
这些方程组一般无法用简单的代数方法求解,需要借助计算机等工具进行求解。
本文将介绍一些常见的多元变量方程组的求解方法。
一、高斯-约旦消元法高斯-约旦消元法是求解线性方程组的一种常见方法,其基本思想是通过多次消元,使方程组限制的范围不断缩小,最终求得方程组的解。
具体步骤如下:1.将方程组写成增广矩阵的形式;2.选定一个系数矩阵的元素作为主元,通常选择第一行第一列元素,即A[1][1];3.对于其他行的该列元素,减去主元所在行对应元素的倍数,使其变为0;4.重复2-3步骤,直到将矩阵化为上三角矩阵;5.从最后一行开始,依次计算出未知变量的值。
高斯-约旦消元法的复杂度为O(n^3),当方程组的规模较大时,求解速度会非常慢。
二、雅可比迭代法雅可比迭代法是通过迭代求解变量的值,直到收敛于方程组的解的方法,其基本思想是将方程组的每个变量下一次迭代时的值,视为其它变量的当前值,通过逐步迭代,求解出未知变量的值。
具体步骤如下:1.将方程组表示为矩阵形式:Ax=b;2.选择一个初值向量x0,设x^(k)为第k次迭代的结果;3.根据迭代公式x_i^(k+1)=[b_i-(sum(A_ij*x_j^(k)))/(A_ii)]/A_ii,计算x^(k+1),其中i表示第i个未知变量,j表示其它未知变量;4.重复3步骤,直到收敛于方程组的解。
雅可比迭代法适用于系数矩阵为对角占优矩阵的情况,当矩阵的条件数较大时,迭代次数可能会非常多,计算速度较慢。
三、列主元高斯消元法列主元高斯消元法是对高斯-约旦消元法的改进,其主要思想是在每次消元时,选择系数矩阵中绝对值最大的元素作为主元,以此来避免出现数值精度过低等问题。
具体步骤如下:1.将方程组写成增广矩阵的形式;2.选定一个未知数作为主元,使得该列元素的绝对值最大;3.将该列中主元所在行交换到最上面;4.对于其他行的该列元素,减去主元所在行对应元素的倍数,使其变为0;5.重复2-3-4步骤,直到将矩阵化为上三角矩阵;6.从最后一行开始,依次计算出未知变量的值。
Guass列选主元消去法和三角分解法

Guass列选主元消去法和三⾓分解法 最近数值计算学了Guass列主消元法和三⾓分解法解线性⽅程组,具体原理如下:1、Guass列选主元消去法对于AX =B1)、消元过程:将(A|B)进⾏变换为,其中是上三⾓矩阵。
即:k从1到n-1a、列选主元选取第k列中绝对值最⼤元素作为主元。
b、换⾏c、归⼀化d、消元2)、回代过程:由解出。
2、三⾓分解法(Doolittle分解)将A分解为如下形式由矩阵乘法原理a、计算U的第⼀⾏,再计算L的第⼀列b、设已求出U的1⾄r-1⾏,L的1⾄r-1列。
先计算U的第r⾏,再计算L的第r列。
a)计算U的r⾏b)计算L的r列C#代码: 代码说明:Guass列主消元法部分将计算出来的根仍然储存在增⼴矩阵的最后⼀列,⽽Doolittle分解,将分解后的结果也储存⾄原来的数组中,这样可以节约空间。
using System;using System.Windows.Forms;namespace Test{public partial class Form1 : Form{public Form1(){InitializeComponent();}private void Cannel_Button_Click(object sender, EventArgs e){this.textBox1.Clear();this.textBox2.Clear();this.textBox3.Clear();boBox1.SelectedIndex = -1;}public double[,] GetNum(string str, int n){string[] strnum = str.Split(' ');double[,] a = new double[n, n + 1];int k = 0;for (int i = 0; i < n; i++){for (int j = 0; j < strnum.Length / n; j++){a[i, j] = double.Parse((strnum[k]).ToString());k++;}}return a;}public void Gauss(double[,] a, int n){int i, j;SelectColE(a, n);for (i = n - 1; i >= 0; i--){for (j = i + 1; j < n; j++)a[i, n] -= a[i, j] * a[j, n];a[i, n] /= a[i, i];}}//选择列主元并进⾏消元public void SelectColE(double[,] a, int n){int i, j, k, maxRowE;double temp; //⽤于记录消元时的因数for (j = 0; j < n; j++){maxRowE = j;for (i = j; i < n; i++)if (System.Math.Abs(a[i, j]) > System.Math.Abs(a[maxRowE, j]))maxRowE = i;if (maxRowE != j)swapRow(a, j, maxRowE, n); //与最⼤主元所在⾏交换//消元for (i = j + 1; i < n; i++){temp = a[i, j] / a[j, j];for (k = j; k < n + 1; k++)a[i, k] -= a[j, k] * temp;}}return;}public void swapRow(double[,] a, int m, int maxRowE, int n){int k;double temp;for (k = m; k < n + 1; k++){temp = a[m, k];a[m, k] = a[maxRowE, k];a[maxRowE, k] = temp;}}public void Doolittle(double[,] a, int n){for (int i = 0; i < n; i++){if (i == 0){for (int j = i + 1; j < n; j++)a[j, 0] = a[j, 0] / a[0, 0];}else{double temp = 0, s = 0;for (int j = i; j < n; j++){for (int k = 0; k < i; k++){temp = temp + a[i, k] * a[k, j];}a[i, j] = a[i, j] - temp;}for (int j = i + 1; j < n; j++){for (int k = 0; k < i; k++){s = s + a[j, k] * a[k, i];}a[j, i] = (a[j, i] - s) / a[i, i];}}}}private void Exit_Button_Click(object sender, EventArgs e){this.Close();}private void Confirm_Button_Click(object sender, EventArgs e){if (this.textBox2.Text.Trim().ToString().Length == 0){this.textBox2.Text = this.textBox1.Text.Trim();}else{this.textBox2.Text = this.textBox2.Text + "\r\n" + this.textBox1.Text.Trim();}this.textBox1.Clear();}private void Calculate_Button_Click(object sender, EventArgs e){string str = this.textBox2.Text.Trim().ToString();string myString = str.Replace("\n", " ").Replace("\r", string.Empty);double[,] a = new double[this.textBox2.Lines.GetUpperBound(0) + 1, this.textBox2.Lines.GetUpperBound(0) + 2];a = GetNum(myString, this.textBox2.Lines.GetUpperBound(0) + 1);if (boBox1.Text == "Guass列主消元法"){Gauss(a, this.textBox2.Lines.GetUpperBound(0) + 1);for (int i = 0; i < this.textBox2.Lines.GetUpperBound(0) + 1; i++){this.textBox3.Text = this.textBox3.Text + "\r\nX" + (i + 1) + "=" + a[i, this.textBox2.Lines.GetUpperBound(0) + 1]; }}else if (boBox1.Text == "Doolittle三⾓分解法"){this.textBox3.Enabled = true;Doolittle(a, this.textBox2.Lines.GetUpperBound(0) + 1);bel3.Text = "分解后的结果:";this.textBox3.Clear();this.textBox3.Text += "L矩阵:\r\n";for (int i = 0; i < this.textBox2.Lines.GetUpperBound(0) + 1; i++) {for (int j = 0; j < this.textBox2.Lines.GetUpperBound(0) + 1; j++) {if (j < i){this.textBox3.Text += a[i, j].ToString() + "\t";}else if (i == j){this.textBox3.Text += "1\t";}else{this.textBox3.Text += "0\t";}}this.textBox3.Text += "\r\n";}this.textBox3.Text += "\r\nU矩阵:\r\n";for (int i = 0; i < this.textBox2.Lines.GetUpperBound(0) + 1; i++) {for (int j = 0; j < this.textBox2.Lines.GetUpperBound(0) + 1; j++) {if (j >= i){this.textBox3.Text += a[i, j].ToString() + "\t";}else{this.textBox3.Text += "0\t";}}this.textBox3.Text += "\r\n";}}}private void textBox1_KeyDown(object sender, KeyEventArgs e){if (e.KeyCode == Keys.Enter){if (this.textBox1.Text.Trim().ToString().Length == 0){Calculate_Button_Click(sender, e);}else{Confirm_Button_Click(sender, e);}}}private void button1_Click(object sender, EventArgs e){this.textBox2.Enabled = true;}}} 运⾏截图: ⾄此完毕。
数值分析--解线性方程组的直接方法

值 为A的特征值,x为A对应的特征向量,A的全体特征值
分 析
称为A的谱,计作 ( A),即 ( A) {i ,i 1,2,, n}, 则称
》
( A)
max
1in
|
i
|
为矩阵A的谱 半 径.
三、特殊矩阵
第5章 解线性方程组的直接方法
1) 对角矩阵
2) 三对角矩阵
3) 上三角矩阵
4) 上海森伯(Hessenberg)阵
分 析
1.00x 1.00y 2.00
》 解法1: 1.00105 x 1.00 y 1.00
(1.00 1.00105) y (2.00 1.00105)
1.00105 x 1.00 y 1.00
1.00
105
y
1.00
105
x 0.00,
y 1.00
第5章 解线性方程组的直接方法
1
Ly b y 3,Ux y x 1.
2
1
第5章 解线性方程组的直接方法
§3 高斯主元素消去法
若ak(kk) 0,或ak(kk)很接近于0,会导致其他元素数量级严重 增长和舍入误差的扩散,使得计算结果不可靠.
《例3’采用3位十进制,用消元法求解
数 值
1.00105 x 1.00y 1.00
L21L1 U2U11
L21L1
U
U 1
21
I
(因为上式右边为上三角矩阵,左边为单位下三角矩阵
从而上式两边都必须等于单位矩阵)
《 数
L1 L2 , U1 U2
1 1 1
值分例2
析
.例1中,A
0
4
-1,将A作LU分解。
范数-摆脱课本繁琐的公式,比较好懂

p
范数的特殊情况。 注:前三种范数都是p—范数的特殊情况。其中 前三种范数都是 范数的特殊情况
|| X ||∞ = lim || X || p
p →∞
计算方法三 计算方法三⑤
向量范数的连续性: 向量范数的连续性
5/35
定理3.3 设f(X)=||X||为Rn上的任一向量范数 则f(X) 定理 为 上的任一向量范数,则 的分量x 的连续函数. 为X的分量 1,x2,…,xn的连续函数 的分量
lim x i = xi (i = 1,2,..., n)
(k ) k →∞
则称向量X= (x1,x2,...,xn)T为向量序列 则称向量 , {X(k)}的极限,或者说向量序列 (k)}收敛 的极限, 的极限 或者说向量序列{X 收敛 于向量X, 于向量 ,记为
lim X
k →∞
(k )
=X 或 X
(k )
→ X (k → ∞)
计算方法三 计算方法三⑤
计算方法三 计算方法三⑤
x1 (k ) ( k ) x2 X = ………… M x (k ) n (k ) x1 x1 (k ) x2 ( k ) x2 X = → = M M x (k ) x n n
几种常用的矩阵范数: 几种常用的矩阵范数:
n
13/35
a11 a21 设 A= ⋅⋅⋅ a n1
a12 ⋅⋅⋅ a1n A 1 = max∑aij 列范数 1≤j≤n i=1 n a22 ⋅⋅⋅ a2n A ∞ = max∑aij 行范数 ⋅⋅⋅ ⋅⋅⋅ ⋅⋅⋅ 1≤i≤n j=1 T an2 ⋅⋅⋅ ann A 2 = λ (A A) max AF =
线性代数求解方法和技巧

线性代数求解方法和技巧线性代数是数学中重要的一个分支,研究向量空间、线性变换和线性方程组等内容。
在实际问题中,我们常常需要用线性代数的方法来解决问题,因此掌握线性代数的求解方法和技巧对于理解和应用数学是非常重要的。
首先,我们讨论线性方程组的求解方法。
线性方程组是由一组线性方程组成的方程组,其中每个方程的未知数的次数都为1。
对于n个未知数和m个方程的线性方程组,我们有以下几种常用的求解方法:1. 列主元消元法:这是最常用的线性方程组求解方法之一。
它的基本思想是通过行变换将线性方程组化为一个三角形式,进而求解得到方程组的解。
在进行行变换时,要选择合适的列主元,即选择主元元素绝对值最大的一列作为主元素。
2. 矩阵求逆法:对于一个可逆的n阶方阵A,我们可以通过求A的逆矩阵来求解线性方程组Ax=b。
具体地,我们首先通过高斯消元法将方程组化为三角形式,然后根据三角形式的矩阵求逆公式来求解x。
3. LU分解法:对于一个n阶非奇异矩阵A,我们可以将其分解为一个下三角矩阵L和一个上三角矩阵U的乘积,即A=LU。
接着,我们可以通过LU分解来求解线性方程组Ax=b。
具体地,我们首先通过LU分解将方程组化为Lc=b和Ux=c两个方程组,然后依次求解这两个方程组得到x的值。
除了以上的求解方法,还有一些线性方程组的特殊情况和对应的求解方法:1. 齐次线性方程组:如果线性方程组右边的常数项都为0,即b=0,那么我们称为齐次线性方程组。
对于齐次线性方程组,其解空间是一个向量空间。
我们可以通过高斯消元法来求解齐次线性方程组,先将其化为三角形式,然后确定自由未知量的个数,最后确定解空间的基底。
2. 奇异线性方程组:如果线性方程组的系数矩阵A是奇异矩阵,即det(A)=0,那么我们称为奇异线性方程组。
对于奇异线性方程组,其解可能不存在,或者存在无穷多解。
我们可以通过计算矩阵A的秩来确定线性方程组的解的情况。
另外,在实际问题中,我们可能会遇到大规模的线性方程组,这时候求解方法和技巧还需要考虑到计算效率的问题。
用列主元消去法解线性方程组

用列主元消去法解线性方程组
列主元消去法是一种解决线性方程组的有效方法,它可以有效地解决多元一次
方程组的求解问题。
列主元消去法是一种基于矩阵分解的算法,它将线性方程组转换为一系列的三
角形矩阵,从而解决线性方程组的求解问题。
它的基本思想是,首先将矩阵A分解为两个矩阵,即列主元矩阵P和三角矩阵U,然后将原方程组PA=LU转换为LUx=Pb,其中Pb是P矩阵乘以b向量的结果,最后求解LUx=Pb即可得到x向量,从而解决线性方程组的求解问题。
列主元消去法的优点在于它可以有效地解决多元一次方程组的求解问题,而且
它的计算量较小,可以在较短的时间内完成求解。
此外,它还可以有效地处理矩阵的奇异性,即当矩阵A的行列式为零时,它仍然可以有效地求解线性方程组。
列主元消去法在互联网领域的应用也很广泛,它可以用于解决复杂的数据挖掘
问题,如推荐系统、搜索引擎等,这些问题都可以用列主元消去法来解决。
此外,它还可以用于解决机器学习中的优化问题,如支持向量机、神经网络等,这些问题也可以用列主元消去法来解决。
总之,列主元消去法是一种有效的解决线性方程组的方法,它可以有效地解决
多元一次方程组的求解问题,而且它的计算量较小,可以在较短的时间内完成求解。
此外,它在互联网领域的应用也很广泛,可以用于解决复杂的数据挖掘问题和机器学习中的优化问题。
实验三 高斯消去法和三角分解法1

实验报告实验三 高斯消去法与矩阵的三角分解一、实验目的1、掌握列主元素消去法,并且能够用MATLAB 编写相关程序,实现高斯消去法的求解。
2、能够用矩阵理论理解与研究高斯消去法,通过对矩阵的初等变换实现高斯消去法。
3、学会矩阵的三角分解,并且能够用MATLAB 编写相关程序,实现矩阵的三角分解,解方程组。
二、上机内容⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡2822171310871234567112345611123451111234111112311111121111111764321x x x x x x1、用列主元素高斯消去法求解方程组。
2、用列主元消去法求解方程组(实现PA=LU) 要求输出: (1)计算解X;(2)L,U;(3)正整型数组IP(i),(i=1,···,n) (记录主行信息)。
三、实验原理1、列主元素消去法用高斯消去法求解方程组时,为了减小误差,在消去的过程中要避免用绝对值较小的主元素。
因此在高斯消去法的每一步应该在系数矩阵货消去后的低阶矩阵中选取绝对值较大的元素作为主元素,保持|m ik |<=1,以减小计算过程中的舍入误差对计算解的影响。
此方法为完全主元素消去法。
完全主元素消去法在选主元素时花费一定的计算机时间,因此实际计算中常用列主元消去法。
列主元消去法在每次选主元时,仅依次按列选取绝对值最大的元素作为主元素,且仅交换两行,再进行消元计算。
装订 线第k步计算如下:对于k=1,2,…,n-1(1)按列选主元:即确定t使(2)如果t≠k,则交换[A,b]第t行与第k行元素。
(3)消元计算(4)回代求解计算流程图回代求解 b=b/a (当a nn ≠0)b ←(b -∑a x )/adet=a nn *det输出计算解及行列式及detk=1,2,…,n-1输入n ,A,b,εdet=1按列主元|a i(k),k |=max|a ik |C 0=a i(k),k换行 a ik a i(k)j(j=k,…n ) b k b j(k), 消元计算 (i=k+1,…,n ) a ik=a ik -a kk *m ik a ij=a ij -a kj *m ik (j=k+1,…,n )|C 0|<εi k =kdet=a kk det否否是是k<=n-1输出det(A)=0停机停机2. 矩阵的三角分解法 (1)定理设 n n R A ⨯∈ 。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
计算方法实验报告1【课题名称】用列主元高斯消去法和列主元三角分解法解线性方程【目的和意义】高斯消去法是一个古老的求解线性方程组的方法,但由它改进得到的选主元的高斯消去法则是目前计算机上常用的解低阶稠密矩阵方程组的有效方法。
用高斯消去法解线性方程组的基本思想时用矩阵行的初等变换将系数矩阵A 约化为具有简单形式的矩阵(上三角矩阵、单位矩阵等),而三角形方程组则可以直接回带求解 用高斯消去法解线性方程组b Ax =(其中A ∈Rn ×n )的计算量为:乘除法运算步骤为32(1)(1)(21)(1)(1)262233n n n n n n n n n n nMD n ----+=+++=+-,加减运算步骤为(1)(21)(1)(1)(1)(25)6226n n n n n n n n n n AS -----+=++=。
相比之下,传统的克莱姆法则则较为繁琐,如求解20阶线性方程组,克莱姆法则大约要19510⨯次乘法,而用高斯消去法只需要3060次乘除法。
在高斯消去法运算的过程中,如果出现abs(A(i,i))等于零或过小的情况,则会导致矩阵元素数量级严重增长和舍入误差的扩散,使得最后的计算结果不可靠,所以目前计算机上常用的解低阶稠密矩阵方程的快速有效的方法时列主元高斯消去法,从而使计算结果更加精确。
2、列主元三角分解法高斯消去法的消去过程,实质上是将A 分解为两个三角矩阵的乘积A=LU ,并求解Ly=b 的过程。
回带过程就是求解上三角方程组Ux=y 。
所以在实际的运算中,矩阵L 和U 可以直接计算出,而不需要任何中间步骤,从而在计算过程中将高斯消去法的步骤进行了进一步的简略,大大提高了运算速度,这就是三角分解法 采用选主元的方式与列主元高斯消去法一样,也是为了避免除数过小,从而保证了计算的精确度【计算公式】1、 列主元高斯消去法设有线性方程组Ax=b ,其中设A 为非奇异矩阵。
方程组的增广矩阵为第1步(k=1):首先在A 的第一列中选取绝对值最大的元素1l a ,作为第一步的主元素:111211212222112[,]n n n l n nn n a a a a b a a a b a a a b ⎡⎤⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦a b然后交换(A ,b )的第1行与第l 行元素,再进行消元计算。
设列主元素消去法已经完成第1步到第k -1步的按列选主元,交换两行,消元计算得到与原方程组等价的方程组 A(k)x=b(k)第k 步计算如下:对于k=1,2,…,n -1(1)按列选主元:即确定t 使 (2)如果t ≠k ,则交换[A ,b]第t 行与第k 行元素。
(3)消元计算消元乘数mik 满足:(4)回代求解2、 列主元三角分解法 对方程组的增广矩阵 经过k -1步分解后,可变成如下形式:111max 0l i i n a a ≤≤=≠(1)(1)(1)(1)(1)1112111(2)(2)(2)(2)22222()(()1)()()()()()1,1()(,)()[,][,] k k k k nk k nk n k k k k k kk kn k k k k n k k k n nn a a a a b a a a b a a b a b b a a a +++⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥→=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦A b A b ()()max 0k k tk ik k i na a ≤≤=≠,(1,,)ik ik ik kka a m i k n a ←=-=+, (,1,,), (1,,)ij ij ik kji i ik k a a m a i j k n b b m b i k n ←+=+⎧⎨←+=+⎩⎪⎪⎩⎪⎪⎨⎧--=-←←∑+=)1,,2,1(,)(1n n i a x a b x a b x ii n i j j ij i i nnn n [,]A A b =11121,11111222,122221,11,1,1,211,11,2121,112,112,1k k k k k k k j n k k j n k k k i i i k n n kk kj kn k ik ij in i nknjk k k j k n n nnk k n a a a b A a u u u u u u y l l l l l l ll l l l u u u u u y u u u u y a a b a a b l a -------------⎡→⎣⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎦第k 步分解,为了避免用绝对值很小的数kku 作除数,引进量1111 (,1,,;1,2,,) ()/ (1,2,,;1,2,,)k kj kj km mj m k ik ik im mk kkm u a l u j k k n k n l a l u u i k k n k n -=-=⎧=-=+=⎪⎪⎨⎪=-=++=⎪⎩∑∑11(,1,,)k i ik im mk m s a l u i k k n -==-=+∑,于是有kk u =ks 。
如果 ,则将矩阵的第t 行与第k 行元素互换,将(i ,j )位置的新元素仍记为jjl 或jja ,然后再做第k 步分解,这时【列主元高斯消去法程序流程图】max t i k i ns s ≤≤= ()/ 1,2,,)1 (1,2,,),kk k k t iki k ik u s s s l s s i k k n l i k k n ===++≤=++即交换前的,(且【列主元高斯消去法Matlab主程序】function x=gauss1(A,b,c) %列主元法高斯消去法解线性方程Ax=bif (length(A)~=length(b)) %判断输入的方程组是否有误disp('输入方程有误!')return;enddisp('原方程为AX=b:') %显示方程组Abdisp('------------------------')n=length(A);for k=1:n-1 %找列主元[p,q]=max(abs(A(k:n,k))); %找出第k列中的最大值,其下标为[p,q]q=q+k-1; %q在A(k:n,k)中的行号转换为在A中的行号if abs(p)<cdisp('列元素太小,det(A)≈0');break;elseif q>ktemp1=A(k,:); %列主元所在行不是当前行,将当前行与列主A(k,:)=A(q,:); 元所在行交换(包括b)A(q,:)=temp1;temp2=b(k,:);b(k,:)=b(q,:);b(q,:)=temp2;end%消元for i=k+1:nm(i,k)=A(i,k)/A(k,k); %A(k,k)将A(i,k)消为0所乘系数A(i,k:n)=A(i,k:n)-m(i,k)*A(k,k:n); %第i行消元处理b(i)=b(i)-m(i,k)*b(k); %b消元处理endenddisp('消元后所得到的上三角阵是')A %显示消元后的系数矩阵b(n)=b(n)/A(n,n); %回代求解for i=n-1:-1:1b(i)=(b(i)-sum(A(i,i+1:n)*b(i+1:n)))/A(i,i);endclear x;disp('AX=b的解x是') x=b;【调用函数解题】【列主元三角分解法程序流程图】【列主元三角分解法Matlab主程序】①自己编的程序:function x=PLU(A,b,eps) %定义函数列主元三角分解法函数if (length(A)~=length(b)) %判断输入的方程组是否有误disp('输入方程有误!')return;enddisp('原方程为AX=b:') %显示方程组Abdisp('------------------------')n=length(A);A=[A b]; %将A与b合并,得到增广矩阵for r=1:nif r==1for i=1:n[c d]=max(abs(A(:,1))); %选取最大列向量,并做行交换if c<=eps %最大值小于e,主元太小,程序结束break;elseendd=d+1-1;p=A(1,:);A(1,:)=A(d,:);A(d,:)=p;A(1,i)=A(1,i);endA(1,2:n)=A(1,2:n);A(2:n,1)=A(2:n,1)/A(1,1); %求u(1,i)elseu(r,r)=A(r,r)-A(r,1:r-1)*A(1:r-1,r); %按照方程求取u(r,i)if abs(u(r,r))<=eps %如果u(r,r)小于e,则交换行p=A(r,:);A(r,:)=A(r+1,:);A(r+1,:)=p;elseendfor i=r:nA(r,i)=A(r,i)-A(r,1:r-1)*A(1:r-1,i); %根据公式求解,并把结果存在矩阵A中endfor i=r+1:nA(i,r)=(A(i,r)-A(i,1:r-1)*A(1:r-1,r))/A(r,r); %根据公式求解,并把结果存在矩阵A中endendendy(1)=A(1,n+1);for i=2:nh=0;for k=1:i-1h=h+A(i,k)*y(k);endy(i)=A(i,n+1)-h; %根据公式求解y(i)endx(n)=y(n)/A(n,n);for i=n-1:-1:1h=0;for k=i+1:nh=h+A(i,k)*x(k);endx(i)=(y(i)-h)/A(i,i); %根据公式求解x(i)endAdisp('AX=b的解x是')x=x'; %输出方程的解②可直接得到P,L,U并解出方程解的的程序(查阅资料得子函数PLU1,其作用是将矩阵A分解成L乘以U的形式。
PLU2为调用PLU1解题的程序,是自己编的)(Ⅰ).function [l,u,p]=PLU1(A) %定义子函数,其功能为列主元三角分解系数矩阵A [m,n]=size(A); %判断系数矩阵是否为方阵if m~=nerror('矩阵不是方阵')returnendif det(A)==0 %判断系数矩阵能否被三角分解error('矩阵不能被三角分解')endu=A;p=eye(m);l=eye(m); %将系数矩阵三角分解,分别求出P,L,Ufor i=1:mfor j=i:mt(j)=u(j,i);for k=1:i-1t(j)=t(j)-u(j,k)*u(k,i);endenda=i;b=abs(t(i));for j=i+1:mif b<abs(t(j))b=abs(t(j));a=j;endendif a~=ifor j=1:mc=u(i,j);u(i,j)=u(a,j);u(a,j)=c;endfor j=1:mc=p(i,j);p(i,j)=p(a,j);p(a,j)=c;endc=t(a);t(a)=t(i);t(i)=c;endu(i,i)=t(i);for j=i+1:mu(j,i)=t(j)/t(i);endfor j=i+1:mfor k=1:i-1u(i,j)=u(i,j)-u(i,k)*u(k,j);endendendl=tril(u,-1)+eye(m);u=triu(u,0)(Ⅱ).function x=PLU2(A,b) %定义列主元三角分解法的函数[l,u,p]=PLU1(A) %调用PLU分解系数矩阵A m=length(A); %由于A左乘p,故b也要左乘p v=b;for q=1:mb(q)=sum(p(q,1:m)*v(1:m,1));endb(1)=b(1) %求解方程Ly=b for i=2:1:mb(i)=(b(i)-sum(l(i,1:i-1)*b(1:i-1)));endb(m)=b(m)/u(m,m); %求解方程Ux=y for i=m-1:-1:1b(i)=(b(i)-sum(u(i,i+1:m)*b(i+1:m)))/u(i,i);endclear x;disp('AX=b的解x是') x=b;【调用函数解题】①②【编程疑难】这是第一次用matlab编程,对matlab的语句还不是非常熟悉,因此在编程过程中,出现了许多错误提示。