Gauss列主元消去法
Gauss列主元素消去法实验
Lab06.Gauss 列主元素消去法实验【实验目的和要求】1.使学生深入理解并掌握Gauss 消去法和Gauss 列主元素消去法步骤; 2.通过对Gauss 消去法和Gauss 列主元素消去法的程序设计,以提高学生程序设计的能力;3.对具体问题,分别用Gauss 消去法和Gauss 列主元素消去法求解。
通过对结果的分析比较,使学生感受Gauss 列主元素消去法优点。
【实验内容】1.根据Matlab 语言特点,描述Gauss 消去法和Gauss 列主元素消去法步骤。
2.编写用不选主元的直接三角分解法解线性方程组Ax=b 的M 文件。
要求输出Ax=b 中矩阵A 及向量b ,A=LU 分解的L 与U ,det A 及解向量x 。
3.编写用Gauss 列主元素消去法解线性方程组Ax=b 的M 文件。
要求输出Ax=b 中矩阵A 及向量b 、PA=LU 分解的L 与U 、det A 及解向量x ,交换顺序。
4.给定方程组(1) ⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡--11134.981.4987.023.116.427.199.103.601.3321x x x(2) ⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡----15900001.582012151********.23107104321x x x x 先用编写的程序计算,再将(1)中的系数3.01改为3.00,0.987改为0.990;将(2)中的系数2.099999改为2.1,5.900001改为9.5,再用Gauss 列主元素消去法解,并将两次计算的结果进行比较。
【实验仪器与软件】1.CPU 主频在1GHz 以上,内存在128Mb 以上的PC ;2.Matlab 6.0及以上版本。
实验讲评:实验成绩:评阅教师:200 年 月 日Lab06.Gauss 列主元素消去法实验第一题:1、算法描述:Ⅰ、Gauss 消去法由书上定理5可知 设Ax=b ,其中A ∈R^(n(1)如果()0(1,2,....,1)k kka k n ≠=-,则可通过高斯消去法将Ax=b 约化为等价的 角形线性方程组,且计算公式为:① 消元计算(k=1,2,….,n-1)()()(1)()()(1)()()/,1,...,,,,1,...,,,1,...,.k k ik ik kk k k k ij ij ik kj k k k iiik k m a a i k n a a m a i j k n b b m b i k n ++==+=-=+=-=+② 回带公式()()()()()1/,()/,1,...,2,1.n n n n nn ni i i i iii j ii j i x b a x ba x a i n =+==-=-∑(2)如果A 为非奇异矩阵,则可通过高斯消去法将方程组Ax=b 约化方程组为上三角矩阵以上消元和回代过程总的乘除法次数为332333nn nn +-≈,加减法次数为32353263nnn n+-≈以上过程就叫高斯消去法。
Gauss列主元消去法、QR(MATLAB)
Gauss列主元消去法、QR(MATLAB)Gauss列主元消去法是一种线性方程组的求解方法,也称Gauss消去法。
其基本思想是将方程组转化为上三角矩阵,然后通过反向代入求解。
该方法的优点在于计算精度高,求解速度快,但缺点是需要大量的计算,尤其是在矩阵阶数较高时。
具体来讲,Gauss列主元消去法的步骤如下:步骤一:将系数矩阵A进行LU分解,其中L是下三角矩阵、U是上三角矩阵。
设$A=LU$,则原方程组可以写成$LUx=b$。
步骤二:通过初等矩阵左乘系数矩阵A,将每一列的主元变为该列所有元素中绝对值最大的那个元素。
这个过程称为选主元,可以避免计算中的数值不稳定问题。
步骤三:将选主元后的系数矩阵A进行LU分解,得到$L^{'}$、$U^{'}$。
步骤五:通过反向代入求解$U^{'}x=y$,得到$x$的解。
Gauss列主元消去法的实现通常通过矩阵的变换来实现。
对于$n$阶矩阵$A=[a_{ij}]$,通过一系列的行变换,可以将其变为上三角矩阵。
其中的变换可以表示为:$$ R_{i} \leftrightarrow R_{j} $$其中,$R_{i}$和$R_{j}$分别表示矩阵$A$中的第$i$行和第$j$行,$k$是一个非零常数。
这些变换被称为初等行变换。
在MATLAB中,可以使用已经实现好的{\color{blue}\texttt{gauss}}函数来求解线性方程组。
该函数实现的算法是Gauss列主元消去法。
其调用格式为:x = gauss(A,b)其中,$A$是系数矩阵,$b$是结果向量。
函数返回结果向量$x$。
如果$A$或$b$不合法,则函数会返回一个空向量。
除了Gauss列主元消去法,还有一种常用的求解线性方程组的方法是QR分解法。
步骤二:通过正交矩阵左乘系数矩阵$A$,使其变为一个上三角矩阵。
这个过程称为正交相似变换。
步骤三:将$b$进行正交相似变换,得到$Q^{T}b$。
高斯列主元消去法
高斯列主元消去法2.3高斯列主元消去法解线性方程组一:问题的提出我们都知道,高斯列主元素消去法是计算机上常用来求解线性方程组的一种直接的方法。
就是在不考虑舍入误差的情况下,经过有限步的四则运算可以得到线性方程组的准确解的一类方法。
实际运算的时候因为只能有限小数去计算,因此只能得到近似值。
在实际运算的时候,我们很多时候也常用高斯消去法。
但是高斯消去法在计算机中运算的时候常会碰到两个问题。
1.一旦遇到一些主元等于0,消元过程便无法进行下去。
2.在长期使用中还发现,即使消元过程能进行下去,但是当一些主元的绝对值很小时,求解出的结果与真实结果相差甚远。
为了避免高斯消去法消元过程中出现的上述两个问题,一般采用所谓的选择主元法。
其中又可以分为列选主元和全面选主元两种方法。
目前计算机上常用的按列选主元的方法。
因此我在这里做的也是列选主元高斯消去法。
二、算法的基本思想大家知道,如果一个线性方程组的系数矩阵是上三角矩阵时,即这种方程组我们称之为上三角方程组,它是很容易求解的。
我们只要把方程组的最下面的一个方程求解出来,在把求得的解带入倒数第二个方程,求出第二个解,依次往上回代求解。
然而,现实中大多数线性方程组都不是上面所说的上三角方程组,所以我们有可以把不是上三角的方程通过一定的算法化成上三角方程组,由此我们可以很方便地求出方程组的解。
高斯消元法的目的就是把一般线性方程组简化成上三角方程组。
于是高斯消元法的基本思想是:通过逐次消元将所给的线性方程组化为上三角形方程组,继而通过回代过程求解线性方程组。
三、算法的描述1、设有n元线性方程组如下:=2、第一步:如果a11!=0,令li1= ai1/a11, I= 2,3,……,n用(-li1)乘第一个方程加到第i个方程上,得同解方程组:a(1)11a(1)12...a(1)1nx1b(1)1a(1)21a(1)22...a(1)2nx2b(1)2 .......=.a(1)n-11 a(1)n-12 . .a(1)n-1nxn-1b(1)n-1a(1)n1a(1)n2. . . a(1)nnxnb(1)n简记为:A(2)x=b(2)其中a(2)ij = a(1)ij – li1 a(1)1j ,I ,j = 2,3,..,nb(2)I = b(1)I – li1 b(1)1 ,I = 2,3,...,n第二步:如果a(2)22!=0,令li2= a(2)i2/a(2)22, I= 3,……,n依据同样的原理,对矩阵进行化间(省略),依次下去,直到完成!最后,得到上三角方程组:a(1)11a(1)12...a(1)1nx1b(1)1a(1)22 ...a(1)2nx2b(1)2 .......=.. .a(n-1)n-1nxn-1b(n-1)n-1. . . a(n)nnxnb(n)n简记为:A(n)x=b(n)最后从方程组的最后一个方程进行回代求解为:n = b(n) / a(n)nni = ( b(k)k - a(k)kjxj ) / a(k)kk以上为高斯消去法的基本过程。
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;}}} 运⾏截图: ⾄此完毕。
2-2 Gauss列主元消去法
S2 若ann 0,则输出“ A是奇异矩阵”;停机 . S3 置xn an,n1 / ann ;
对i n 1, n 2,...1,
ai,n1 n aij x j
置xi
j i 1
aii
S4 输出x1, x2,..., xn ;停机.
作业:
P50 习题3
k in
aik
;
S12 若aik ,k 0,则输出“ A是奇异矩阵”;停机 .
S13 若ik k,则
akj aik , j j k,...,n 1;
S14 对i k 1,..., n
置aik aik / akk ; 对j k 1,..., n 1
置aij aij aik akj.
§2-2 Gauss列主元消去法
一、Gauss列主元消去法的引入 例1. 用3位浮点数运算,求解线性方程组
0.0001xx11
x2 x2
1 2
解: 本方程组的精度较高的解为
x* (1.00010001 ,0.99989999 )T
用Gauss消去法求解
A ( A,b)
0.000100 1
1 1
21
0.000100
m2110 000
0
回代后得到
1
1
1.00 104 1.00 104
x1 0.00 , x2 1.00
与精确解相比,该结果显然是错误的 究其原因,在求行乘数时用了很小的数0.0001作除数
如果在求解时将1,2行交换,即
A ( A,b)
1 0.000100
1 1
a(2) i2
,
交换第2行和第i2行,
2in
然后进行消元,得[ A(3) , b(3) ].
列主元素Gauss消去法Jacobi迭代法原理及计算方法
一、 列主元素Gauss 消去法、Jacobi 迭代法原理及计算方法1. 列主元素Gauss 消去法:1.1 Gauss 消去法基本原理设有方程组Ax b =,设A 是可逆矩阵。
高斯消去法的基本思想就是将矩阵的初等行变换作用于方程组的增广矩阵[]B A b = ,将其中的A 变换成一个上三角矩阵,然后求解这个三角形方程组。
1.2 列主元Gauss 消去法计算步骤将方程组用增广矩阵[]()(1)ijn n B A b a ⨯+== 表示。
1). 消元过程对1,2,,1k n =-(1) 选主元,找{},1,,k i k k n ∈+ 使得 ,max k i k ik k i na a ≤≤= (2) 如果,0k i k a =,则矩阵A 奇异,程序结束;否则执行(3)。
(3) 如果k i k ≠,则交换第k 行与第k i 行对应元素位置,k kj i j a a ↔,,,1j k n =+ 。
(4) 消元,对,,i k n = ,计算/,ik ik kk l a a =对1,,1j k n =++ ,计算.ij ij ik kj a a l a =-2). 回代过程(1) 若0,nn a =则矩阵奇异,程序结束;否则执行(2)。
(2) ,1/;n n n nn x a a +=对1,,2,1i n =- ,计算,11/n i i n ij j ii j i x a a x a +=+⎛⎫=- ⎪⎝⎭∑2. Jacobi 迭代法2.1 Jacobi 迭代法基本原理Jacobi 迭代法的基本思想是对n 元线性方程组b Ax =,.,n n R b R A ∈∈将其变形为等价方程组f Bx x +=,其中.,,n n n n R x R f R B ∈∈∈⨯B 成为迭代矩阵。
从某一取定的初始向量)0(x 出发,按照一个适当的迭代公式 ,逐次计算出向量f Bx x k k +=+)()1( ( 1,0=k ),使得向量序列}{)(k x 收敛于方程组的精确解.(1)输入1,,,,)0(=k n xb A ε,. (2) )(1,1)0()1(∑≠=-=n j i i j ij i iii x a b a x )1,0(n i = (3)判断 ε≤--≤≤)0()1(10max i i n i x x ,若是,输出1)1(2)1(1,,n x x x ,若否,置1+=k k ,)1()0(i i x x =,)2,1(n i =。
Gauss列主元消去法程序设计
《Gauss列主元消去法》实验报告实验名称:Gauss列主元消去法程序设计成绩:___________专业班级:数学与应用数学1202班姓名:王晓阳学号:28实验日期:2014 年11月10日实验报告日期:2014年11月10日一.实验目的1.学习Gauss消去法的基本思路和迭代步骤.2.学会运用matlab编写高斯消去法和列主元消去法程序,求解线性方程组.3.当()ka绝对值较小时,采用高斯列主元消去法.kk4.培养编程与上机调试能力.二、实验内容用消去法解线性方程组的基本思想是用逐次消去未知数的方法把原线性方程组Ax b=化为与其等价的三角形线性方程组,而求解三角形线性方程组可用回代的方法求解.1.求解一般线性方程组的高斯消去法.(1)消元过程:设()0k kk a ≠,第i 个方程减去第k 个方程的()()/k k ik ik kk m a a =倍,(1,,)i k n =+L ,得到()()11k k A x b ++=.()()()()()()()11,,1,,k k k ij ij ik kj k k k ii ik k a a m a i j k n b b m b ++⎧=-=+⎪⎨=-⎪⎩L 经过n-1次消元,可把方程组()()11A x b =化为上三角方程组()()n n A x b =.(2)回代过程:()()()()()1//,1,,1n n n n nn n i i i i i ij j ii j i x b a x b a x a i n =+⎧=⎪⎛⎫⎨=-=- ⎪⎪⎝⎭⎩∑L 以解如下线性方程组为例测试结果.1212312310773264556x x x x x x x x -=⎧⎪-++=⎨⎪-+=⎩2.列主元消去法由高斯消去法可知,在消元过程中可能出现()0k kk a =的情况,这是消去法将无法进行,即使主元素()0k kk a ≠但很小时,用其作除数,会导致其他元素数量级的严重增长和舍入误差的扩散,最后也使得计算解不可靠.这时就需要选取主元素,假定线性方程组的系数矩阵A 是菲奇异的.(1)消元过程:对于1,2,,1k n =-L ,进行如下步骤:1) 按列选主元,记max pk ik k i na a ≤≤=2) 交换增广阵A 的p,k 两行的元素。
课件:5.2(b) Gauss列主元消去法
L L11L21 Ln11为单位下三角矩阵 U A(n)为上三角矩阵
即 且
1
1,1 ln,1
ln1,2 ln,2
ln1,2 ln,3
1 ln,n1
1
a(1) 11
U A(n)
a(1) 12
a(2) 22
a(1) 1n
a(2) 2n
a(n) nn
1 2 3
解: 这个方程组和例1一样,若用Gauss消去法计算会有 小数作除数的现象,若采用换行的技巧,则可避免
108
A (A,b) 1
2
2
3 1
3.712 4.623 2
1.072
5.64
3
3
108 很小,绝对值最大
的列元素为a13 2, 因此1,3行交换
r1 r3
2 1
例2所用的方法是在Gauss消去法的基础上,利用换行 避免小主元作除数,该方法称为Gauss列主元消去法
二、Gauss消元过程与系数矩阵的分解 1.Gauss消去法消元过程的矩阵描述
( A(1) ,b(1) )
a(1) 11
a(1) 21
a(1) 12
a(1) 22
a(1) 1n
a(1) 2n
主元
A (A,b)
0.000100 1
1 1
1 2
l21 10000
0.000100 0
1 1.00 104
1 1.00 104
回代后得到
x1 0.00 , x2 1.00
与精确解相比,该结果相当糟糕 究其原因,在求乘数时用了很小的数0.0001作除数
如果在求解时将1,2行交换,即
A (A,b)
det Ak 0 k 1,2,, n
gauss列主元素消去法matlab
高斯列主元素消去法是一种解线性方程组的常用方法,特别在数值分析和线性代数中应用广泛。
在Matlab中,我们可以使用该方法来解决大规模的线性方程组,包括矩阵的求解和矩阵的反转。
一、高斯列主元素消去法的基本原理高斯列主元素消去法是一种基于矩阵消元的方法,它通过一系列的矩阵变换将原始的线性方程组转化为上三角形式,然后再进行回代求解。
这个方法的核心就是通过矩阵的变换来简化原始的线性方程组,使得求解过程更加简单高效。
在Matlab中,我们可以利用矩阵运算和函数来实现高斯列主元素消去法,如`lu`分解函数和`\"`运算符等。
通过这些工具,我们能够快速地求解各种规模的线性方程组并得到准确的结果。
二、高斯列主元素消去法在Matlab中的实现在Matlab中,我们可以通过调用`lu`函数来实现高斯列主元素消去法。
该函数返回一个上三角矩阵U和一个置换矩阵P,使得PA=LU。
通过对U进行回代求解,我们可以得到线性方程组的解。
除了`lu`函数之外,Matlab还提供了一些其他的函数和工具来帮助我们实现高斯列主元素消去法,比如`\"`运算符和`inv`函数等。
通过这些工具的组合使用,我们能够更加灵活地进行线性方程组的求解,并且可以方便地处理特殊情况和边界条件。
三、高斯列主元素消去法的应用与局限性高斯列主元素消去法在实际应用中具有广泛的适用性,特别是对于大规模的线性方程组或者稀疏矩阵的求解。
通过Matlab中的工具和函数,我们可以快速地求解各种规模的线性方程组,并得到高精度的数值解。
然而,高斯列主元素消去法也存在一些局限性,比如对于奇异矩阵或者接近奇异矩阵的情况时,该方法的求解精度可能会下降。
在实际应用中,我们需要结合具体的问题和矩阵特性来选择合适的求解方法,以确保得到准确的结果。
四、个人观点和总结作为一种经典的线性方程组求解方法,高斯列主元素消去法在Matlab 中具有较好的实现和应用效果。
通过对其原理和实现细节的深入理解,我们能够更加灵活地应用该方法,并且能够更好地理解其适用性和局限性。
Gauss列主元消去法
贵州师范大学数学与计算机科学学院学生实验报告课程名称: 数值分析 班级: 实验日期: 年 月 日 学 号: 姓名: 指导教师: 实验成绩:一、实验名称实验五:线性方程组的数值解法二、实验目的及要求1. 让学生掌握用列主元gauss 消去法、超松弛迭代法求解线性方程组.2. 培养Matlab 编程与上机调试能力.三、实验环境每人一台计算机,要求安装Windows XP 操作系统,Microsoft office2003、MATLAB6.5(或7.0).四、实验内容1. 编制逐次超松弛迭代(SOR 迭代)函数(子程序),并用于求解方程组⎪⎪⎩⎪⎪⎨⎧=-++=+-+=++-=+++-141414144321432143214321x x x x x x x x x x x x x x x x取初始向量T x )1,1,1,1()0(=,迭代控制条件为 5)1()(1021||||--⨯≤-k k x x 请绘制出迭代次数与松弛因子关系的函数曲线,给出最佳松弛因子.SOR 迭代的收敛速度是否一定比Gauss-Seidel 迭代快?2. 编制列主元 Gauss 消去法函数(子程序),并用于解 ⎪⎩⎪⎨⎧=++-=-+-=+-615318153312321321321x x x x x x x x x要求输出方程组的解和消元后的增广矩阵. 注:题2必须写实验报告五、算法描述及实验步骤Gauss 消去法:功能 解方程组b Ax = .输入 n ,n n ij a A ⨯=)(,T n b b b b ),,,(21 =.输出 方程组的解T n x x x x ),,,(21 =或失败信息.步1 对1,,2,1-=n k 执行步2→步4 .步2 调选列主元模块 .步3 若0=kk a ,则=x “消去法失败”,结束 .步4 对n k k i ,,2,1 ++=执行步5→步6 .步5 对n k k j ,,2,1 ++=执行ij kj kk ik ij a a a a a +⨯-⇐/ .步6 i k kk ik i b b a a b +⨯-⇐/ .步7 nn n n a b x /⇐ .步8 对1,,2,1 --=n n i 执行ii n i j j ij i i a x a b x /)(1∑+=-⇐ .步9 输出T n x x x x ),,,(21 = .选列主元模块:功能 选列主元 .输入 n k k i b n k k j i a i ij ,,1,,;,,1,,, +=+= .输出 n k k i b n k k j i a i ij ,,1,,;,,1,,, +=+= .步1 kk a m ⇐;k l ⇐ .步2 对n k k i ,,2,1 ++=执行若m a ik >则ik a m ⇐;i l ⇐ .步3 若k l ≠,则交换kj a 和lj a ,n k k j ,,1, +=;交换k b 和l b .步4 返回主模块 .六、调试过程及实验结果>> A=[12,-3,3;-18,3,-1;1,1,1];>> b=[15;-15;6];>> x=Gauss1(A,b)Ab =-18.0000 3.0000 -1.0000 -15.00000 1.1667 0.9444 5.16670 0 3.1429 9.4286 index = 1x = 1.0000 2.0000 3.0000七、总结由于数)1(-k kka 在Gauss 消去法中有着突出的作用,第k 步消元时,要用)1(-k kk a 作除数,如果)1(-k kk a =0消元会失败,即使主元)1(-k kk a ≠0,但很小时,舍入误差也会使计算结果面目全非,避免这种缺陷的基本方法就是选主元。
实验二:Gauss列主元消去法
实验二:Gauss列主元消去法程序1:Gauss列主元消去法A=input('请输入线性方程组的增广矩阵A=');n=length(A)-1;x=zeros(n,1);aa=zeros(n,1);for j=1:nfor i=1:(n+1)AA(j,i)=abs(A(j,i));endendfor k=1:(n-1)for i=k:naa(i-(k-1))=AA(i,k);endfor i=k:nif AA(i,k)==max(aa)breakendendif AA(i,k)==0breakfprintf('方程组系数矩阵奇异\n');elsefor j=k:(n+1)jh=A(i,j);A(i,j)=A(k,j);A(k,j)=jh;endendfenzi=A(k,k);for j=k:(n+1)A(k,j)=A(k,j)/fenzi;endfor p=(k+1):njj=A(p,k);for j=k:(n+1)A(p,j)=A(p,j)-jj*A(k,j);endendendif k==(n-1)x(n)=A(n,(n+1))/A(n,n);for i=(n-1):(-1):1he=0;for j=(i+1):nhe=he+A(i,j)*x(j);endx(i)=A(i,(n+1))-he;endendx用Gauss列主元消去法解方程组:1.请输入线性方程组的增广矩阵A=[1e-008,2,3,1;-1,3.172,4.623,2;-2,1.072,5. 643,3]x =-0.4653-0.07000.38002.请输入线性方程组的增广矩阵A=[4,-2,4,10;-2,17,10,3;-4,10,9,-7];x =2.94640.6071-0.14293.请输入线性方程组的增广矩阵A=[0.3e-020,1,0.7;1,1,0.9]x =0.20000.7000程序2:不选主元的高斯消去法A=input('请输入线性方程组的增广矩阵A=');n=length(A)-1;x=zeros(n,1);for k=1:(n-1)if A(k,k)==0breakfprintf('方程组不能用普通的高斯消去法解\n');elsefenzi=A(k,k);for j=k:(n+1)A(k,j)=A(k,j)/fenzi;endfor p=(k+1):njj=A(p,k);for j=k:(n+1)A(p,j)=A(p,j)-jj*A(k,j);endendx(n)=A(n,(n+1))/A(n,n);for i=(n-1):(-1):1he=0;for j=(i+1):nhe=he+A(i,j)*x(j);endx(i)=A(i,(n+1))-he;endendendx用不选主元的Gauss消去法解方程组:1.请输入线性方程组的增广矩阵A=[4,-2,4,10;-2,17,10,3;-4,10,9,-7];x =2.94640.6071-0.14292.请输入线性方程组的增广矩阵A=[1e-008,2,3,1;-1,3.172,4.623,2;-2,1.072,5. 643,3];x =-0.4653-0.07000.38003.请输入线性方程组的增广矩阵A=[0.3e-020,1,0.7;1,1,0.9]x =0.7000。
数学实验题目5 相对Gauss列主元消去法
数学实验题目5 相对Gauss 列主元消去法摘要由一般线性方程组在使用Gauss 消去法求解时,从求解过程中可以清楚地看到,若(1)0k kk a -=,必须施以行交换的手续,才能使消去过程继续下去。
有时既使(1)0k kk a -≠,但其绝对值很小,由于舍入误差的影响,消去过程也会出现不稳定现象。
因此,为使这种不稳定现象发生的可能性减至最小,在施行消去过程时每一步都要选主元素,即要寻找行r ,使(1)(1)||max ||k k rk ik i ka a -->=并将第r 行与第k 行交换,以使(1)k kka -的当前值(即(1)k ika -的数值)远大于0。
这种列主元消去法的主要步骤如下: 1.消元过程 对1,2,,1k n =-,做1º 选主元,记||max ||rk ik i ka a >=若0rk a =,说明方程组系数矩阵奇异,则停止计算,否则进行2º。
2º 交换A (增广矩阵)的,r k 两行元素,,1rj kja a j k n ↔=+3º 计算/ij ij ik kj kk a a a a a =-1,,i k n =+1,,1j k n =++2.回代过程 对,1,,2,1k n n =-,计算,11(/)nk k n kjj kk j k x a ax a +=-=-∑前言利用Gauss 列主元消去法、显式相对Gauss 列主元消去法、隐式相对Gauss 列主元消去法求解线性方程组程序设计流程是否否是开 始输入A (增广矩阵)1k = ||max ||rk ik i ka a >=rk a =1,,1,,1/ij ij ik kj kk i k n j k n a a a a a =+=++=-交换A 中,r k 两行1k n <-,11,1,,2,1()/nk k n kjj kkj k k n n x a ax a +=+=-=-∑输出x结 束 1k =问题1(1)程序运行如下:x = GaussSysSolve(Mat1_1,b1_1)x = 1.0000 1.0000 1.0000 1.0000x = GaussExpSysSolve(Mat1_1,b1_1)x = 1.0000 1.0000 1.0000 1.0000x = GaussIneSysSolve(Mat1_1,b1_1)x = 1.0000 1.0000 1.0000 1.0000 (2)程序运行如下:x = GaussSysSolve(Mat1_2,b1_2)x = 1.0000 1.0000 1.0000 1.0000x = GaussExpSysSolve(Mat1_2,b1_2)x = 1.0000 1.0000 1.0000 1.0000x = GaussIneSysSolve(Mat1_2,b1_2)x = 1.0000 1.0000 1.0000 1.0000 (3)程序运行如下:x = GaussSysSolve(Mat1_3,b1_3)x = 1.0000 1.0000 1.0000 1.0000x = GaussExpSysSolve(Mat1_3,b1_3)x = 1.0000 1.0000 1.0000 1.0000x = GaussIneSysSolve(Mat1_3,b1_3)x = 1.0000 1.0000 1.0000 1.0000 (4)程序运行如下:x = GaussSysSolve(Mat1_4,b1_4)x = 1.0000 1.0000 1.0000 1.0000x = GaussExpSysSolve(Mat1_4,b1_4)x = 1.0000 1.0000 1.0000 1.0000x = GaussIneSysSolve(Mat1_4,b1_4)x = 1.0000 1.0000 1.0000 1.0000问题2(1)程序运行如下:= GaussSysSolve(Mat2_1,b2_1)x = 1.0915 0.2832 1.1463 -0.1008x = GaussExpSysSolve(Mat2_1,b2_1)x = 1.0915 0.2832 1.1463 -0.1008x = GaussIneSysSolve(Mat2_1,b2_1)x = 1.0915 0.2832 1.1463 -0.1008 (2)程序运行如下:x = GaussSysSolve(Mat2_2,b2_2)x = 0.5162 0.4152 0.1100 1.0365x = GaussExpSysSolve(Mat2_2,b2_2)x = 0.5162 0.4152 0.1100 1.0365x = GaussIneSysSolve(Mat2_2,b2_2)x = 0.5162 0.4152 0.1100 1.0365 (3)程序运行如下:x = GaussSysSolve(Mat2_3,b2_3)x = 1.0000 1.0000 1.0000x = GaussExpSysSolve(Mat2_3,b2_3)x = 1 1 1x = GaussIneSysSolve(Mat2_3,b2_3)x = 1.0000 1.0000 1.0000(4)程序运行如下:x = GaussSysSolve(Mat2_4,b2_4)x = 1 1 1x = GaussExpSysSolve(Mat2_4,b2_4)x = 1.0000 1.0000 1.0000x = GaussIneSysSolve(Mat2_4,b2_4)x = 1 1 1使用的函数function x = GaussSysSolve(A, b)% GaussSysSolve 用Gauss消去法解线性方程组Ax = b%% Synopsis: x = GaussSysSolve(A, b)%% Input: A = 系数矩阵% b = 方程组右端%% Output: x = 线性系统的解向量[m,n] = size(A);b = b(:); %将b变为列向量if m ~= n %A必须为方阵error('Argument matrix A must be square!');elseif m ~= length(b) %b的长度应与A维度相同error('The dimentions of A and b do not agree!');endAb = [A b]; %构造增广矩阵for i = 1:n[amax, imax] = max(Ab(i:n, i)); %选择主元if amax == 0 %主元为0,矩阵奇异error('Tne Linear System is singular!');elseif i ~= imax+i-1 %主元行数与i不同时,交换这两行Ab([i imax+i-1],:) = Ab([imax+i-1 i], :);endfor j = i+1:n %向下消元Ab(j,:) = Ab(j,:) - Ab(i,:) * Ab(j,i)/amax;endendx = zeros(n,1);x(n) = Ab(n,n+1)/Ab(n,n);for k = n-1:-1:1 %计算x x(k) = ( Ab(k,n+1) - Ab(k,k+1:n)*x(k+1:n) ) / Ab(k,k);endfunction x = GaussExpSysSolve(A, b)% GaussExpSysSolve 用显式Gauss列主元消去法解线性方程组Ax = b%% Synopsis: x = GaussExpSysSolve(A, b)%% Input: A = 系数矩阵% b = 方程组右端%% Output: x = 线性系统的解向量[m,n] = size(A);b = b(:); %将b变为列向量if m ~= n %A必须为方阵error('Argument matrix A must be square!');elseif m ~= length(b) %b的长度应与A维度相同error('The dimentions of A and b do not agree!');endAb = [A b]; %构造增广矩阵for i = 1:n %显式平衡技术s = max(Ab(i,1:n));Ab(i,:) = Ab(i,:)/s;endfor i = 1:n[amax, imax] = max(Ab(i:n, i)); %选择主元if amax == 0 %主元为0,矩阵奇异error('Tne Linear System is singular!');elseif i ~= imax+i-1 %主元行数与i不同时,交换这两行Ab([i imax+i-1],:) = Ab([imax+i-1 i], :);endfor j = i+1:n %向下消元Ab(j,:) = Ab(j,:) - Ab(i,:) * Ab(j,i)/amax;endendx = zeros(n,1);x(n) = Ab(n,n+1)/Ab(n,n);for k = n-1:-1:1 %计算xx(k) = ( Ab(k,n+1) - Ab(k,k+1:n)*x(k+1:n) ) / Ab(k,k);endfunction [x,det] = GaussIneSysSolve(A, b)% GaussIneSysSolve 用隐式Gauss列主元消去法解线性方程组Ax = b%% Synopsis: x = GaussIneSysSolve(A, b)%% Input: A = 系数矩阵% b = 方程组右端%% Output: x = 线性系统的解向量% det = 系数矩阵行列式的值[m,n] = size(A);b = b(:); %将b变为列向量if m ~= n %A必须为方阵error('Argument matrix A must be square!');elseif m ~= length(b) %b的长度应与A维度相同error('The dimentions of A and b do not agree!');endAb = [A b]; %构造增广矩阵det = 1; %初始化系数矩阵行列式为1for i = 1:n %隐式平衡技术s(i) = max(abs(Ab(i,1:n)));if s(i) == 0error('Tne Linear System is singular!'); %系数矩阵某行全为0时,矩阵奇异endends = s(:);for k = 1:n-1[c, kmax] = max(abs(Ab(k:n, k)./s(k:n))); %选择主元if c == 0 %主元为0,矩阵奇异det = 0;error('Tne Linear System is singular! det(A) = 0');elseif k ~= kmax+k-1 %主元行数与k不同时,交换这两行s([k kmax+k-1]) = s([kmax+k-1 k]);Ab([k kmax+k-1],:) = Ab([kmax+k-1 k], :);det = -det;endfor j = k+1:n %向下消元Ab(j,:) = Ab(j,:) - Ab(k,:) * Ab(j,k)/Ab(k,k);enddet = Ab(k,k)*det;endif Ab(n,n) == 0det = 0;error('Tne Linear System is singular!'); %最后一行唯一非0元素为0时,矩阵奇异endx = zeros(n,1);x(n) = Ab(n,n+1)/Ab(n,n);for k = n-1:-1:1 %计算x x(k) = ( Ab(k,n+1) - Ab(k,k+1:n)*x(k+1:n) ) / Ab(k,k);enddet = Ab(n,n)*det;思考题(1)在各主元不是非常小的时候,三种方法结果一致(2)隐式平衡列选主元法最好,应为当主元很小时,普通的Gauss消元法会产生很大误差;显式平衡列选主元法每一行除以其绝对值最大元素时会引入额外的舍入误差。
解线性方程组的列主元素高斯消去法和lu分解法
数值试验报告分析一、实验名称:解线性方程组的列主元素高斯消去法和LU 分解法二、实验目的及要求:通过数值实验,从中体会解线性方程组选主元的必要性和LU分解法的优点,以及方程组系数矩阵和右端向量的微小变化对解向量的影响。
三、算法描述:本次试验采用的是高斯列主元消去法和LU分解法求解线性方程组的解。
其中,高斯消去法的基本思想是避免接近于零的数作分母;能进行到底的条件: 当A可逆时,列主元Gauss(高斯)消去法一定能进行到底。
优点: 具有很好的数值稳定性;具有与顺序Gauss 消去法相同的计算量。
列主元Gauss(高斯)消去法的精度显著高于顺序Gauss(高斯)消去法。
注意:省去换列的步骤,每次仅选一列中最大的元。
矩阵的三角分解法是A=LU,L 是下三角阵,U是上三角阵,Doolittle 分解:L 是单位下三角阵,U是上三角阵;Crout 分解:L 是下三角阵,U是单位上三角阵。
矩阵三角分解的条件是矩阵 A 有唯一的Doolittle 分解的充要条件是 A 的前n-1 顺序主子式非零;矩阵A有唯一的Crout 分解的充要条件是 A 的前n-1 顺序主子式非零。
三角分解的实现是通过(1)Doolittle 分解的实现;(2)Doolittle 分解的缺点:条件苛刻,且不具有数值稳定性。
(3)用Doolittle 分解求解方程组: AX=b LUX=b LY=bA=LU UX=Y ;四、实验内容:解下列两个线性方程组3.01 6.03 1.99 x1 11) 1.27 4.16 1.23 x2 10.987 4.81 9.34 x3 110 7 0 1 x1 83 2.099999 6 2 x2 5.9000012) 5 1 5 1x3 52 1 0 2 x4 1a、用你熟悉的算法语言编写程序用列主元高斯消去法和LU分解求解上述两个方程组,输出Ax=b 中矩阵 A 及向量b, A=LU 分解的L 及U,detA 及解向量x.b、将方程组(1)中系数 3.01 改为 3.00 ,0.987 改为0.990 ,用列主元高斯消去法求解变换后的方程组,输出列主元行交换次序,解向量x 及detA ,并与(1)中结果比较。
Matlab程序Gauss列主元消去法
Matlab程序Gauss列主元消去法编辑整理:尊敬的读者朋友们:这里是精品文档编辑中心,本文档内容是由我和我的同事精心编辑整理后发布的,发布之前我们对文中内容进行仔细校对,但是难免会有疏漏的地方,但是任然希望(Matlab程序Gauss列主元消去法)的内容能够给您的工作和学习带来便利。
同时也真诚的希望收到您的建议和反馈,这将是我们进步的源泉,前进的动力。
本文可编辑可修改,如果觉得对您有帮助请收藏以便随时查阅,最后祝您生活愉快业绩进步,以下为Matlab程序Gauss列主元消去法的全部内容。
.Gauss 列主元消去法求解线性方程组迭代法计算停止的条件为:.采用用用Gauss 列主元消去法时,Matlab 计算程序为: clearclcA=[2 2 1 2;4 1 3 —1;—4 —2 0 1;2 3 2 3];B=[1;2;1;0];n=length(B);X=zeros (n,1);c=zeros(1,n);d1=0;for i=1:n —1max=abs(A(i ,i ));m=i;for j=i+1:nif max<abs (A (j,i ))max=abs (A (j,i));m=j;endendif (m ~=i)for k=i :nc(k)=A(i ,k);A (i,k)=A(m ,k);A (m,k)=c (k );endd1=B(i);12346212425027,2085113270x x x x -⎛⎫⎛⎫⎛⎫ ⎪ ⎪ ⎪- ⎪ ⎪ ⎪= ⎪ ⎪ ⎪-- ⎪ ⎪ ⎪⎝⎭⎝⎭⎝⎭6)()1(3110max -+≤≤<-k j k j j x xB(i)=B(m);B(m)=d1;endfor k=i+1:nfor j=i+1:nA(k,j)=A(k,j)—A(i,j)*A(k,i)/A(i,i);endB(k)=B(k)—B(i)*A(k,i)/A(i,i);A(k,i)=0;endendX(n)=B(n)/A(n,n);for i=n-1:—1:1sum=0;for j=i+1:nsum=sum+A(i,j)*X(j);endX(i)=(B(i)-sum)/A(i,i);endX计算结果为:X =(1。
高斯列主元消去法
问题提出:采用高斯列主元消去法解线性方程组。
算法(公式)推导:高斯顺序消去法有一个最大的缺点就是一旦对角元素为0,就进行不下去了,为了解决这个问题就有了高斯主元消去法。
如果在高斯顺序消去法消去过程进行到第i 步时,先选取a ri ()n r i ≤≤中(即第i 列)绝对值最大的元素,设为第j 行的元素aji ,然后将第i+1行至第n 行中的每一行减去第i 行乘以ii kj a a (k 代表行号),依次进行消元,这样得到的算法叫高斯按列主元消去法。
高斯按列主元消去法的算法步骤介绍如下:1. 将方程组写成以下的增广矩阵的形式: 432144434241343332312423222114131211b b b b a a a a a a a a a a a a a a a a 2. 对k=1,2,3,…..,n-1,令∑==nk s sk pk a a max ,交换增广矩阵的第k 行与第p 行;对j=k+1,K+2,……..,n,计算*km jkjm jm kk a a a a a =-(m=k,k+1,....n)kk jk k j j a a b b b *-=算法结束。
3. 在MATLABE 中编程实现的高斯按列主元消去法函数为:GaussXQLineMain功能:高斯按列主元消去法求线性方程组Ax=b 的解调用格式:[x,XA]=GaussXQLineMain(A,b)其中,A :线性方程组的系数矩阵;B:线性方程组中的常数向量;x:线性方程组的解:XA:消元后的系数矩阵(可选的输出参数)。
高斯列主元消去法用MATLAB实现如下所示:4.其中用到上三角矩阵求解函数:在MATLABE中编程实现的上三角系数矩阵求解函数为:SolveUPTriangle 功能:求上三角系数矩阵的线性方程组Ax=b的解调用格式:x=SolveUpTriangel(A,b)其中,A :线性方程组的系数矩阵;b :线性方程组中的常数向量; X :线性方程组的解;上三角系数矩阵求解函数用MATLAB 实现如下所示:高斯按列主元消去法解线性方程组应用实例:用高斯按列主元消去法求解下列线性方程组的解。
解线性方程组的列主元素高斯消去法和lu分解法
数值试验报告分析一、实验名称:解线性方程组的列主元素高斯消去法和LU 分解法二、实验目的及要求:通过数值实验,从中体会解线性方程组选主元的必要性和LU 分解法的优点,以及方程组系数矩阵和右端向量的微小变化对解向量的影响。
三、算法描述:本次试验采用的是高斯列主元消去法和LU 分解法求解线性方程组的解。
其中,高斯消去法的基本思想是避免接近于零的数作分母;能进行到底的条件:当A 可逆时,列主元Gauss(高斯)消去法一定能进行到底。
优点:具有很好的数值稳定性;具有与顺序Gauss 消去法相同的计算量。
列主元Gauss(高斯)消去法的精度显著高于顺序Gauss(高斯)消去法。
注意:省去换列的步骤,每次仅选一列中最大的元。
矩阵的三角分解法是A=LU,L 是下三角阵,U 是上三角阵,Doolittle 分解:L 是单位下三角阵,U 是上三角阵;Crout 分解:L 是下三角阵,U 是单位上三角阵。
矩阵三角分解的条件 是矩阵A 有唯一的Doolittle 分解的充要条件是A 的前n-1顺序主子式非零;矩阵A 有唯一的Crout 分解的充要条件是A 的前n-1顺序主子式非零。
三角分解的实现是通过(1)Doolittle 分解的实现; (2)Doolittle 分解的缺点:条件苛刻,且不具有数值稳定性。
(3)用Doolittle 分解求解方程组: AX=b LUX=b LY=bA=LU UX=Y ;四、实验内容:解下列两个线性方程组(1) ⎪⎪⎪⎭⎫ ⎝⎛=⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎭⎫ ⎝⎛--11134.981.4987.023.116.427.199.103.601.3321x x x (2) ⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=⎪⎪⎪⎪⎪⎭⎫ ⎝⎛⎪⎪⎪⎪⎪⎭⎫ ⎝⎛----15900001.582012151526099999.23107104321x x x x a 、用你熟悉的算法语言编写程序用列主元高斯消去法和LU 分解求解上述两个方程组,输出Ax=b 中矩阵A 及向量b, A=LU 分解的L 及U ,detA 及解向量x.b 、将方程组(1)中系数3.01改为3.00,0.987改为0.990,用列主元高斯消去法求解变换后的方程组,输出列主元行交换次序,解向量x及detA,并与(1)中结果比较。
数值分析实验作业,gauss消去法的数值稳定性分析
实验3.1 Gauss 消去法的数值稳定性试验实验目的:观察和理解Gauss 消元过程中出现小主元(即)(k kka 很小)时引起的方程组解的数值不稳定性。
实验内容:求解方程组b Ax =,其中(1)⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡⨯=11212592.1121-130.6-291.51314.59103.015-1A ,⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=2178.4617.591b ; (2)⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡----=2010151526990999999999.23107102A ,⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=15019000000000.582b .实验要求:(1) 计算矩阵的条件数,判断系数矩阵是良态的还是病态的。
(2) 用Gauss 列主元消去法求得L 和U 及解向量421,R x x ∈.(3) 用不选主元的Gauss 消去法求得L ~和U ~及解向量421~,~R x x ∈.(4) 观察小主元并分析其对计算结果的影响.程序如下:计算矩阵条件数及Gauss 列主元消去法:format longengA1=[0.3e-15 59.14 3 1;5.291 -6.130 -1 2;11.2 9 5 2;1 2 1 1]; b1=[59.17;46.78;1;2]; n=4;k2=cond(A1) %k2为矩阵的条件数;for k=1:n-1a=max(abs(A1(k:n,k))); [p,k]=find(A1==a); B=A1(k,:);c=b1(k);A1(k,:)=A1(p,:);b1(k)=b1(p); A1(p,:)=B;b1(p)=c; if A1(k,k)~=0A1(k+1:n,k)=A1(k+1:n,k)/A1(k,k);A1(k+1:n,k+1:n)=A1(k+1:n,k+1:n)-A1(k+1:n,k)*A1(k,k+1:n); else break end endL1=tril(A1,0); for i=1:n L1(i,i)=1; end L=L1U=triu(A1,0) for j=1:n-1b1(j)=b1(j)/L(j,j);b1(j+1:n)=b1(j+1:n)-b1(j)*L(j+1:n,j); endb1(n)=b1(n)/L(n,n); for j=n:-1:2b1(j)=b1(j)/U(j,j);b1(1:j-1)=b1(1:j-1)-b1(j)*U(1:j-1,j); endb1(1)=b1(1)/U(1,1); x1=b1运行结果如下: K2=68.43;⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡---⨯=-14929.00202.00893.0011755.04724.00011079.2600118L⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡-=801.0000231.1835.2001314.5902592.11U 1x =[18.9882;3.3378;-34.747;-33.9865] 不选主元的Gauss 消去法程序:clearformat longengA1=[0.3e-15 59.14 3 1;5.291 -6.130 -1 2;11.2 9 5 2;1 2 1 1]; b1=[59.17;46.78;1;2]; n=4;for k=1:n-1A1(k+1:n,k)=A1(k+1:n,k)/A1(k,k);A1(k+1:n,k+1:n)=A1(k+1:n,k+1:n)-A1(k+1:n,k)*A1(k,k+1:n); endL1=tril(A1,0); for i=1:n L1(i,i)=1; end L=L1U=triu(A1,0) for j=1:n-1b1(j)=b1(j)/L(j,j);b1(j+1:n)=b1(j+1:n)-b1(j)*L(j+1:n,j); endb1(n)=b1(n)/L(n,n); for j=n:-1:2b1(j)=b1(j)/U(j,j);b1(1:j-1)=b1(1:j-1)-b1(j)*U(1:j-1,j); endb1(1)=b1(1)/U(1,1); x1=b1程序运行结果如下:⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡⨯⨯⨯=10189.010333.3011168.21033.3700110637.170001~151515L⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡--⨯-⨯-⨯-⨯=-5.000816010637.171091.5210043.101314.59103.0~15151815U ]0;0;0005.1;6848.23[~1=x同理可得2A 对应的系数矩阵条件数及Gauss 列主元消去法求解结果: K2=8.994;⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡-⨯=1333.04.0001104.0-3.0-0015.0000112-L ⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡--=36667.300030.26005.155.2010710U ]0.1;;0.1;0.1;10444.0[152-⨯=-x不选主元的Gauss 消去法结果:⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡⨯⨯=1400.0109999.0-001104998.2-5.00013.0-0001~1212L ⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡⨯⨯⨯--=-3667.3000107495.5109987.14003.26100.1010710~121212U ]000145.1;99994.0;000.1;1045.1[~52-⨯-=-x实验4.5 三次样条插值函数的收敛性问题提出:多项式插值不一定收敛的,即插值的节点多,效果不一定就好。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
贵州师范大学数学与计算机科学学院学生实验报告
课程名称: 数值分析 班级: 数本(一)班 实验日期: 年 月 日
学 号: 0098(81) 姓名: 吴胜 指导教师: 杨一都
实验成绩:
一、实验名称
实验五:线性方程组的数值解法
二、实验目的及要求
1. 让学生掌握用列主元gauss 消去法、超松弛迭代法求解线性方程组.
2. 培养Matlab 编程与上机调试能力.
三、实验环境
每人一台计算机,要求安装Windows XP 操作系统,Microsoft office2003、(或.
四、实验内容
1. 编制逐次超松弛迭代(SOR 迭代)函数(子程序),并用于求解方程组
⎪⎪⎩⎪⎪⎨⎧=-++=+-+=++-=+++-1
41
414144321432143214321x x x x x x x x x x x x x x x x
取初始向量T x )1,1,1,1()0(=,迭代控制条件为 5)1()(102
1||||--⨯≤
-k k x x 请绘制出迭代次数与松弛因子关系的函数曲线,给出最佳松弛因子.SOR 迭代
的收敛速度是否一定比Gauss-Seidel 迭代快
2. 编制列主元 Gauss 消去法函数(子程序),并用于解 ⎪⎩⎪⎨⎧=++-=-+-=+-615318153312321
321321x x x x x x x x x
要求输出方程组的解和消元后的增广矩阵. 注:题2必须写实验报告
五、算法描述及实验步骤
Gauss 消去法:
功能 解方程组b Ax = .
输入 n ,n n ij a A ⨯=)(,T n b b b b ),,,(21 =.
输出 方程组的解T n x x x x ),,,(21 =或失败信息.
步1 对1,,2,1-=n k 执行步2→步4 .
步2 调选列主元模块 .
步3 若0=kk a ,则=x “消去法失败”,结束 .
步4 对n k k i ,,2,1 ++=执行步5→步6 .
步5 对n k k j ,,2,1 ++=执行ij kj kk ik ij a a a a a +⨯-⇐/ .
步6 i k kk ik i b b a a b +⨯-⇐/ .
步7 nn n n a b x /⇐ .
步8 对1,,2,1 --=n n i 执行ii n i j j ij i i a x a b x /)(1∑+=-
⇐ .
步9 输出T n x x x x ),,,(21 = .
选列主元模块:
功能 选列主元 .
输入 n k k i b n k k j i a i ij ,,1,,;,,1,,, +=+= .
输出 n k k i b n k k j i a i ij ,,1,,;,,1,,, +=+= .
步1 kk a m ⇐;k l ⇐ .
步2 对n k k i ,,2,1 ++=执行若m a ik >则ik a m ⇐;i l ⇐ .
步3 若k l ≠,则交换kj a 和lj a ,n k k j ,,1, +=;交换k b 和l b .
步4 返回主模块 .
六、调试过程及实验结果
>> A=[12,-3,3;-18,3,-1;1,1,1];
>> b=[15;-15;6];
>> x=Gauss1(A,b)
Ab =
0 0 index = 1
x =
七、总结
由于数)1(-k kk
a 在Gauss 消去法中有着突出的作用,第k 步消元时,要用)1(-k kk a 作除数,如果)1(-k kk a =0消元会失败,即使主元)1(-k kk a ≠0,但很小时,舍入误差也会使计算结果面目全非,避免这种缺陷的基本方法就是选主元。
通过选主元,就可避免绝对值小的数作除数,从而避免舍入误差的恶性增长,使得Gauss 列主元消去法是解中小规模的线性方程组和某些大型稀疏线性方程组的有效方法。
八、附录(源程序清单)
function [x,index]=Gauss1(A,b)
[n,m]=size(A);x=zeros(n,1); index=1
for k=1:n-1
a_max=0;
for i=k:n
if abs(A(i,k))>a_max
a_max=abs(A(i,k));r=i;
end
end
if a_max<1e-10
index=0;return ;
end
if r>k
for j=k:n
z=A(k,j);A(k,j)=A(r,j);A(r,j)=z;
end
z=b(k);b(k)=b(r);b(r)=z;
end
for i=k+1:n
m=A(i,k)/A(k,k);
for j=k:n
A(i,j)=A(i,j)-m*A(k,j);
end
b(i)=b(i)-m*b(k);
end
end
if abs(A(n,n))==0
index=0;return ;
end
Ab=[A,b]
x(n)=b(n)/A(n,n);
for i=n-1:-1:1
for j=i+1:n
b(i)=b(i)-A(i,j)*x(j);
end
x(i)=b(i)/A(i,i);
end。