计算方法--解线性方程组的直接法
文章目录
- 一、Gauss 消元法
- 1.顺序高斯消元法
- 总计算量
- 2.主元素高斯消元法
- 列主元素高斯消元法
- 3.高斯-约当(Gauss-Jordan)消去法
- 总计算量
- 二、矩阵三角分解法
- 1.直接三角分解法(LU分解、Doolittle分解)
- 条件:
- 计算顺序
- 计算公式
- 计算量
- 优点
- 三、Cholesky分解法(平方根法)
- 条件
- 分解结果
- 计算L的步骤
- Step1,计算对角元素
- Step2,计算同列其他元素
- 改进的Cholesky分解法
- Step 1. 计算 D 值
- Step 2. 计算第j列的L值
- 追赶法
- 条件
- 求解过程
- 误差分析
- 性质
一、Gauss 消元法
Gauss消元法较为简单,只会简单叙述。
1.顺序高斯消元法
通过对方程组的增广矩阵进行初等行变换将矩阵变为上三角矩阵。
再通过回代求得原方程组的解。
总计算量
13n3+n2−13n2\frac{1}{3}n^3+n^2-\frac{1}{3}n^231n3+n2−31n2
2.主元素高斯消元法
为了控制舍入误差的扩大和传播而提出的。
列主元素高斯消元法
将该列绝对值尽可能大的系数作为第k步消元的主元素
3.高斯-约当(Gauss-Jordan)消去法
每次消去对角线下方和上方的元素
每次消元过程中,选取列主元素。
总计算量
12n3+12n2\frac{1}{2}n^3+\frac{1}{2}n^221n3+21n2
计算量比高斯消去法计算量大,但不需要回代。
二、矩阵三角分解法
1.直接三角分解法(LU分解、Doolittle分解)
条件:
矩阵A的各阶顺序主子式都不为0(可逆)。
计算顺序
- 先计算第一行,从左到右依次计算 uuu ,yyy。
- 再计算第一列,从上到下计算 lll。
- 按照前两部依次计算第二行到最后一行。
计算公式
设矩阵 AAA 的增广矩阵 AAA 为 A[n+1][n+2];A[n+1][n+2];A[n+1][n+2]; 下标为0的元素位置为0.
A=(a11a12⋯a1n∣b1a21a22⋯a2n∣b2⋮⋮⋮⋮an1an2⋯ann∣bn)\boldsymbol{A}=\left(\begin{array}{cccc} a_{11} & a_{12} & \cdots & a_{1 n} & |& b_1 \\ a_{21} & a_{22} & \cdots & a_{2 n} & |& b_2 \\ \vdots & \vdots & & \vdots & & \vdots \\ a_{n 1} & a_{n 2} & \cdots & a_{n n} & |& b_n \end{array}\right)A=⎝⎜⎜⎜⎛a11a21⋮an1a12a22⋮an2⋯⋯⋯a1na2n⋮ann∣∣∣b1b2⋮bn⎠⎟⎟⎟⎞
则有:
for(int i=1;i<=n;i++){//i表示计算的层数//对于第i层//先计算u,y值,用u值替换a值,用y值替代b值for(int k=i;k<=n+1;k++){for(int r = 1; r < i; r++)A[i][k] = A[i][k] - A[i][r]*A[r][k];}//再计算l值,替代下三角的a值for(int k = i + 1; k <= n; k++){for(int r = 1; r<i; r++){A[k][i] = A[k][i] - A[k][r]*A[r][i];}A[k][i] /= A[k][k];}}
替换后的值如下:
L&U&y=(u11u12u13⋯u1n∣y1l21u22u23⋯u2n∣y2l31l32u33⋯u3n∣y3⋮⋮⋮⋮⋮ln1ln2⋯⋯unn∣yn)\boldsymbol{L\&U\&y}=\left(\begin{array}{ccccc} u_{11} & u_{12} & u_{13} & \cdots & u_{1n} & |&y_1\\ l_{21} & u_{22} & u_{23} & \cdots & u_{2n} & |&y_2\\ l_{31} & l_{32} & u_{33} & \cdots & u_{3n} & |&y_3\\ \vdots & \vdots & \vdots & & \vdots & &\vdots\\ l_{n 1} & l_{n 2} & \cdots & \cdots & u_{nn} & |&y_n \end{array}\right)L&U&y=⎝⎜⎜⎜⎜⎜⎛u11l21l31⋮ln1u12u22l32⋮ln2u13u23u33⋮⋯⋯⋯⋯⋯u1nu2nu3n⋮unn∣∣∣∣y1y2y3⋮yn⎠⎟⎟⎟⎟⎟⎞
及:u34(黑色)−=红色部分的计算值u_{34}(黑色) -= 红色部分的计算值u34(黑色)−=红色部分的计算值(y值同理)
及:
l43(黑色)−=红色部分的计算值l_{43}(黑色) -= 红色部分的计算值l43(黑色)−=红色部分的计算值
然后,
l43(黑色)/=橙色(橙色为对角线上的值)l_{43}(黑色) /= 橙色(橙色为对角线上的值)l43(黑色)/=橙色(橙色为对角线上的值)
最后通过计算Ux=yUx = yUx=y得出xxx。用xxx的值替代yyy值。
有
for(k=i+1;k<n;k++)for(k=i+1;k<n;k++)for(k=i+1;k<n;k++)
A[i][n+1]=A[i][n+1]−A[i][k]∗A[k][n+1]A[i][n+1] = A[i][n+1] - A[i][k]*A[k][n+1]A[i][n+1]=A[i][n+1]−A[i][k]∗A[k][n+1]
A[i][n+1]/=A[i][i]A[i][n+1] /=A[i][i]A[i][n+1]/=A[i][i]
即 :黑色=黑色−红色的积和黑色 = 黑色 - 红色的积和黑色=黑色−红色的积和
再计算:黑色=黑色/橘色黑色 = 黑色 / 橘色黑色=黑色/橘色
计算量
也是 13n3\frac{1}{3}n^331n3数量级,与高斯消元法计算量基本相同。
优点
若是有多个系数矩阵相同而右边项不同的一系列方程组,用直接三角分解法更简单。
三、Cholesky分解法(平方根法)
条件
- 系数矩阵AAA对称,即A=ATA = A^TA=AT
- 系数矩阵AAA正定,即AAA的特征值均为正值。
分解结果
可以将AAA唯一的分解为A=LLTA = LL^TA=LLT
L=(l110⋯0l21l22⋯0⋮⋮⋱⋮ln1ln2⋯lnn)L=\left(\begin{array}{cccc} l_{11} & 0 & \cdots & 0 \\ l_{21} & l_{22} & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ l_{n 1} & l_{n 2} & \cdots & l_{n n} \end{array}\right)L=⎝⎜⎜⎜⎛l11l21⋮ln10l22⋮ln2⋯⋯⋱⋯00⋮lnn⎠⎟⎟⎟⎞
计算L的步骤
从外层向内侧依次计算。同LULULU分解。
Step1,计算对角元素
lkk=akk−∑p=1k−1lkp2\begin{array}{c} l_{k k}=\sqrt{a_{k k}-\sum_{p=1}^{k-1} l_{k p}^{2}} \end{array}lkk=akk−∑p=1k−1lkp2
- 对角元素=对角元素−红色部分的积和对角元素 = 对角元素 - 红色部分的积和对角元素=对角元素−红色部分的积和(注意:A为对称矩阵。即:a14=a41a_{14}=a_{41}a14=a41)
- 对角元素=对角元素对角元素 = \sqrt{对角元素}对角元素=对角元素
Step2,计算同列其他元素
lik=(aik−∑p=1k−1liplkp)/lkk,i=k+1,k+2,⋯,n\begin{array}{c} l_{i k}=\left(a_{i k}-\sum_{p=1}^{k-1} l_{i p} l_{k p}\right) / l_{k k}, \\ i=k+1, k+2, \cdots, n \end{array}lik=(aik−∑p=1k−1liplkp)/lkk,i=k+1,k+2,⋯,n
- A[i][j]=A[i][j]−红色部分的积和A[i][j]= A[i][j]- 红色部分的积和A[i][j]=A[i][j]−红色部分的积和(注意:因为A13并未更新为l13,所以要将A[1][3]替换为A[3][1]这是与UL分解不同的地方因为A_{13}并未更新为l_{13},所以要将A[1][3]替换为A[3][1]这是与UL分解不同的地方因为A13并未更新为l13,所以要将A[1][3]替换为A[3][1]这是与UL分解不同的地方)
- A[i][j]=A[i][j]/对角元素(橘色)A[i][j] = A[i][j]/对角元素(橘色)A[i][j]=A[i][j]/对角元素(橘色)
在完成 Cholesky 分解后, 我们可分别求解以下系数矩阵。即U=LTU = L^TU=LT:
LY=b,LTX=Y\boldsymbol{L} \boldsymbol{Y}=\boldsymbol{b}, \quad \boldsymbol{L}^{\mathrm{T}} \boldsymbol{X}=\boldsymbol{Y}LY=b,LTX=Y
{y1=b1yi=bi−∑k=1i−1likyk,i=2,3,⋯,n\left\{\begin{array}{l} y_{1}=b_{1} \\ y_{i}=b_{i}-\sum_{k=1}^{i-1} l_{i k} y_{k}, \quad i=2,3, \cdots, n \end{array}\right.{y1=b1yi=bi−∑k=1i−1likyk,i=2,3,⋯,n
{xn=yn/unnxi=(yi−∑k=i+1nuikxk)/uii,i=n−1,n−2,⋯,1.\left\{\begin{array}{l} x_{n}=y_{n} / u_{n n} \\ x_{i}=\left(y_{i}-\sum_{k=i+1}^{n} u_{i k} x_{k}\right) / u_{i i}, \quad i=n-1, n-2, \cdots, 1 . \end{array}\right.{xn=yn/unnxi=(yi−∑k=i+1nuikxk)/uii,i=n−1,n−2,⋯,1.
由于上述线性方程组有大量开方运算,故改进Cholesky分解法。
改进的Cholesky分解法
将AAA分解为:
A=LDLTA = LDL^TA=LDLT
Step 1. 计算 D 值
dj=ajj−∑k=1j−1ljkvjk,vjk=ljkdk,d_{j}=a_{j j}-\sum_{k=1}^{j-1} l_{j k} v_{j k}, \quad v_{j k}=l_{j k} d_{k},dj=ajj−k=1∑j−1ljkvjk,vjk=ljkdk,
Step 2. 计算第j列的L值
注意:L为单位下三角矩阵。即:L矩阵的对角元素为1.
lij=(aij−∑k=1j−1likvjk)/dj,i=j+1,j+2,⋯,n.\begin{array}{c} \\ l_{i j}=\left(a_{i j}-\sum_{k=1}^{j-1} l_{i k} v_{j k}\right) / d_{j}, \quad i=j+1, j+2, \cdots, n . \end{array}lij=(aij−∑k=1j−1likvjk)/dj,i=j+1,j+2,⋯,n.
在完成分解后, 我们可分别求解下列系数矩。
LY=b,DLTX=Y\boldsymbol{L Y}=\boldsymbol{b}, \quad \boldsymbol{D} \boldsymbol{L}^{\mathrm{T}} \boldsymbol{X}=\boldsymbol{Y}LY=b,DLTX=Y
追赶法
条件
- 三对角非线性方程组
- 对角占优(∣bi∣>=∣ai∣+∣ci∣|b_i|>=|a_i|+|c_i|∣bi∣>=∣ai∣+∣ci∣)
(b1c1a2b2c2⋱⋱⋱aibici⋱⋱⋱an−1bn−1cn−1anbn)(x1x2⋮xi⋮xn−1xn)=(d1d2⋮di⋮dn−1dn)\left(\begin{array}{ccccccc} b_{1} & c_{1} & & & & & \\ a_{2} & b_{2} & c_{2} & & & & \\ & \ddots & \ddots & \ddots & & & \\ & & a_{i} & b_{i} & c_{i} & & \\ & & & \ddots & \ddots & \ddots & \\ & & & & a_{n-1} & b_{n-1} & c_{n-1} \\ & & & & & a_{n} & b_{n} \end{array}\right)\left(\begin{array}{c} x_{1} \\ x_{2} \\ \vdots \\ x_{i} \\ \vdots \\ x_{n-1} \\ x_{n} \end{array}\right)=\left(\begin{array}{c} d_{1} \\ d_{2} \\ \vdots \\ d_{i} \\ \vdots \\ d_{n-1} \\ d_{n} \end{array}\right)⎝⎜⎜⎜⎜⎜⎜⎜⎜⎛b1a2c1b2⋱c2⋱ai⋱bi⋱ci⋱an−1⋱bn−1ancn−1bn⎠⎟⎟⎟⎟⎟⎟⎟⎟⎞⎝⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎛x1x2⋮xi⋮xn−1xn⎠⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎞=⎝⎜⎜⎜⎜⎜⎜⎜⎜⎜⎜⎛d1d2⋮di⋮dn−1dn⎠⎟⎟⎟⎟⎟⎟⎟⎟⎟⎟⎞
求解过程
- 通过 Gauss 消元法将其化为系数矩阵为单位上三角形矩阵的方程组
(1β11β2⋱⋱1βn−11)(x1x2⋮xn−1xn)=(γ1γ2⋮γn−1γn)\left(\begin{array}{ccccc} 1 & \beta_{1} & & & \\ & 1 & \beta_{2} & & \\ & & \ddots & \ddots & \\ & & & 1 & \beta_{n-1} \\ & & & & 1 \end{array}\right)\left(\begin{array}{c} x_{1} \\ x_{2} \\ \vdots \\ x_{n-1} \\ x_{n} \end{array}\right)=\left(\begin{array}{c} \gamma_{1} \\ \gamma_{2} \\ \vdots \\ \gamma_{n-1} \\ \gamma_{n} \end{array}\right)⎝⎜⎜⎜⎜⎛1β11β2⋱⋱1βn−11⎠⎟⎟⎟⎟⎞⎝⎜⎜⎜⎜⎜⎛x1x2⋮xn−1xn⎠⎟⎟⎟⎟⎟⎞=⎝⎜⎜⎜⎜⎜⎛γ1γ2⋮γn−1γn⎠⎟⎟⎟⎟⎟⎞ - 回代求得方程组的解。
xn=γn,xi=γi−βixi+1,i=n−1,n−2,⋯,2x_{n}=\gamma_{n}, \quad x_{i}=\gamma_{i}-\beta_{i} x_{i+1}, \quad i=n-1, n-2, \cdots, 2xn=γn,xi=γi−βixi+1,i=n−1,n−2,⋯,2
误差分析
cond(A)=∣∣A−1∣∣⋅∣∣A∣∣cond(A) = ||A^{-1}||\cdot||A||cond(A)=∣∣A−1∣∣⋅∣∣A∣∣刻画了方程组Ax=bAx=bAx=b的病态程度,及解对A,bA,bA,b扰动的敏感程度。
性质
对可逆矩阵A,有:
- cond(A)=cond(A−1)cond(A) = cond(A^{-1})cond(A)=cond(A−1)
- cond(kA)=cond(A)cond(kA) = cond(A)cond(kA)=cond(A)
- cond2(U)=1cond_2(U) = 1cond2(U)=1
- 等等
计算方法--解线性方程组的直接法相关推荐
- 【计算方法】解线性方程组的直接法
0.WARNINGS 本文章主要适用于计算方法代码的实现参考,由于本人是python究极小白,为了实验课速成了一些内容,因此会包含较多的暴力解法orz,如有代码错误.可优化的地方欢迎各位大佬指出,感激 ...
- 解线性方程组的直接法
总结自<数值分析>李乃成版 解线性方程组的直接法 高斯消去法 高斯消去法就是初等数学的消元法,这里是将求解步骤规范化,使之便于编写计算机程序(设置好每一步执行的顺序,使程序具有最小的时间复 ...
- matlab解方程实验,MATLAB实验一解线性方程组的直接法
MATLAB实验一解线性方程组的直接法 实 验 报 告 课程名称 数值分析 实验项目 解线性方程组的直接法 专业班级 姓 名 学 号 指导教师 成 绩 日 期 月 日 一. 实验目的 1. 掌握程序的 ...
- 乔利斯基三角分解_解线性方程组的直接法4.1-2.ppt
您所在位置:网站首页 > 海量文档  > 高等教育 > 微积分 解线性方程组的直接法4.1-2.ppt24页 本文档一 ...
- 数值分析 解线性方程组的直接法(一)
数值分析第三章 引言 高斯消去法 高斯消去法的基本思想 n元线性方程组的高随消去法 高斯消去法的计算步骤 列主元高斯消去法 列主元高斯消去法的基本思想 算法步骤 列主元高斯消元法的代码实例 引言 对于 ...
- Python解线性方程组的直接法(6)————求解三对角方程组的追赶法
求解三对角方程组的追赶法 import numpy as npdef zuiganfa(A, d):n = A.shape[0]l = np.mat(np.zeros(n, dtype=float
- 直接法 matlab,解线性方程组直接方法matlab用法.doc
解线性方程组直接方法matlab用法 在这章中我们要学习线性方程组的直接法,特别是适合用数学软件在计算机上求解的方法. 2.1 方程组的逆矩阵解法及其MATLAB程序 2.1.3 线性方程组有解的判定 ...
- matlab lu分解求线性方程组_计算方法(二)直接三角分解法解线性方程组
封面是WH2里春希在编辑部的上司麻理前辈,有一说一,这条线的第一次H有点恶趣味,不是很喜欢. 一:概述 矩阵分解我学过的挺多种,比如极分解,谱分解,满秩分解,正交三角分解还有这里的直接三角分解大部分我 ...
- 计算方法(三)平方根法及其改进解线性方程组
一:概述 本篇文章介绍解线性方程组的平方根法及改进平方根法,适用范围为系数矩阵为正定Hermite矩阵(下称H阵)的线性方程组.这个方法的理论依据我觉得是来自Schur引理的H阵结构定理,从这个角度我 ...
- 计算方法(二)直接三角分解法解线性方程组
一:概述 矩阵分解我学过的挺多种,比如极分解,谱分解,满秩分解,正交三角分解还有这里的直接三角分解大部分我都没有具体运用的经验.但是这里的三角分解的应用就很直白了,就是把矩阵分解为规律的三角矩阵后,我 ...
最新文章
- uniapp实现页面左右滑动,上下滑动事件
- Go基础编程:作用域
- 卷积神经网络原理_怎样设计最优的卷积神经网络架构?| NAS原理剖析
- Firefox beta 开始原生支持 Windows 10 ARM64
- 黄聪:SQL server 2005高可用性之----数据库镜像
- Matlab的部分文件操作
- peewee创建mysql_python – peewee MySQL,如何创建包装SQL构建的ins的自定义字段类型?...
- php复选框样式,如何自定义checkbox样式?附代码
- 2019年2月数据库流行度排行: PostgreSQL攀至历史新高
- 第五十三天:优化网站的常用方法
- 唐诗辑注 —— 辛夷坞、南园十三首、问六十九
- Python基础---时间模块 (二)
- python 除法总返回浮点
- 评人工智能如何走向新阶段?
- keil4如何设置自动缩进_在Keil中 自动格式化 代码
- python selenium 保存网页_使用python/selenium保存完整的网页(包括css、图像)
- JS 中删除节点的两个方法
- android模拟器游戏大全,安卓模拟器游戏大全_小鸡模拟器
- 花菁染料(cas773041-79-5|cas427882-78-8|cas14134-81-7)结构图及合成路线图
- 【5G系列】MICO学习总结(2)