计算方法--解线性方程组的直接法

2023-11-03

一、Gauss 消元法

Gauss消元法较为简单,只会简单叙述。

1.顺序高斯消元法

通过对方程组的增广矩阵进行初等行变换将矩阵变为上三角矩阵
再通过回代求得原方程组的解。

总计算量

1 3 n 3 + n 2 − 1 3 n 2 \frac{1}{3}n^3+n^2-\frac{1}{3}n^2 31n3+n231n2

2.主元素高斯消元法

为了控制舍入误差的扩大和传播而提出的。

列主元素高斯消元法

将该列绝对值尽可能大的系数作为第k步消元的主元素

3.高斯-约当(Gauss-Jordan)消去法

每次消去对角线下方和上方的元素
每次消元过程中,选取列主元素。

总计算量

1 2 n 3 + 1 2 n 2 \frac{1}{2}n^3+\frac{1}{2}n^2 21n3+21n2
计算量比高斯消去法计算量大,但不需要回代。

二、矩阵三角分解法

1.直接三角分解法(LU分解、Doolittle分解)

条件:

矩阵A的各阶顺序主子式都不为0(可逆)。

计算顺序

  • 先计算第一行,从左到右依次计算 u u u y y y
  • 再计算第一列,从上到下计算 l l l
  • 按照前两部依次计算第二行到最后一行。

计算公式

设矩阵 A A A 的增广矩阵 A A A A [ n + 1 ] [ n + 2 ] ; A[n+1][n+2]; A[n+1][n+2]; 下标为0的元素位置为0.
A = ( a 11 a 12 ⋯ a 1 n ∣ b 1 a 21 a 22 ⋯ a 2 n ∣ b 2 ⋮ ⋮ ⋮ ⋮ a n 1 a n 2 ⋯ a n n ∣ b n ) \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=a11a21an1a12a22an2a1na2nannb1b2bn

则有:

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 = ( u 11 u 12 u 13 ⋯ u 1 n ∣ y 1 l 21 u 22 u 23 ⋯ u 2 n ∣ y 2 l 31 l 32 u 33 ⋯ u 3 n ∣ y 3 ⋮ ⋮ ⋮ ⋮ ⋮ l n 1 l n 2 ⋯ ⋯ u n n ∣ y n ) \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=u11l21l31ln1u12u22l32ln2u13u23u33u1nu2nu3nunny1y2y3yn

及: u 34 ( 黑 色 ) − = 红 色 部 分 的 计 算 值 u_{34}(黑色) -= 红色部分的计算值 u34()=(y值同理)
在这里插入图片描述
及:
l 43 ( 黑 色 ) − = 红 色 部 分 的 计 算 值 l_{43}(黑色) -= 红色部分的计算值 l43()=
然后,
l 43 ( 黑 色 ) / = 橙 色 ( 橙 色 为 对 角 线 上 的 值 ) l_{43}(黑色) /= 橙色(橙色为对角线上的值) l43()/=线
在这里插入图片描述

最后通过计算 U x = y Ux = y Ux=y得出 x x x。用 x x x的值替代 y y y值。

f o r ( 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]
即 : 黑 色 = 黑 色 − 红 色 的 积 和 黑色 = 黑色 - 红色的积和 =
再计算: 黑 色 = 黑 色 / 橘 色 黑色 = 黑色 / 橘色 =/
在这里插入图片描述

计算量

也是 1 3 n 3 \frac{1}{3}n^3 31n3数量级,与高斯消元法计算量基本相同。

优点

若是有多个系数矩阵相同而右边项不同的一系列方程组,用直接三角分解法更简单。

三、Cholesky分解法(平方根法)

条件

  • 系数矩阵 A A A对称,即 A = A T A = A^T A=AT
  • 系数矩阵 A A A正定,即 A A A的特征值均为正值。

分解结果

可以将 A A A唯一的分解为 A = L L T A = LL^T A=LLT
L = ( l 11 0 ⋯ 0 l 21 l 22 ⋯ 0 ⋮ ⋮ ⋱ ⋮ l n 1 l n 2 ⋯ l n n ) 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=l11l21ln10l22ln200lnn

计算L的步骤

从外层向内侧依次计算。同 L U LU LU分解。

Step1,计算对角元素

l k k = a k k − ∑ p = 1 k − 1 l k p 2 \begin{array}{c} l_{k k}=\sqrt{a_{k k}-\sum_{p=1}^{k-1} l_{k p}^{2}} \end{array} lkk=akkp=1k1lkp2

  1. 对 角 元 素 = 对 角 元 素 − 红 色 部 分 的 积 和 对角元素 = 对角元素 - 红色部分的积和 =(注意:A为对称矩阵。即: a 14 = a 41 a_{14}=a_{41} a14=a41
  2. 对 角 元 素 = 对 角 元 素 对角元素 = \sqrt{对角元素} =
    在这里插入图片描述

Step2,计算同列其他元素

l i k = ( a i k − ∑ p = 1 k − 1 l i p l k p ) / l k k , 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=(aikp=1k1liplkp)/lkk,i=k+1,k+2,,n

  1. A [ i ] [ j ] = A [ i ] [ j ] − 红 色 部 分 的 积 和 A[i][j]= A[i][j]- 红色部分的积和 A[i][j]=A[i][j](注意: 因 为 A 13 并 未 更 新 为 l 13 , 所 以 要 将 A [ 1 ] [ 3 ] 替 换 为 A [ 3 ] [ 1 ] 这 是 与 U L 分 解 不 同 的 地 方 因为A_{13}并未更新为l_{13},所以要将A[1][3]替换为A[3][1]这是与UL分解不同的地方 A13l13,A[1][3]A[3][1]UL
  2. A [ i ] [ j ] = A [ i ] [ j ] / 对 角 元 素 ( 橘 色 ) A[i][j] = A[i][j]/对角元素(橘色) A[i][j]=A[i][j]/
    在这里插入图片描述

在完成 Cholesky 分解后, 我们可分别求解以下系数矩阵。即 U = L T U = L^T U=LT:

L Y = b , L T X = Y \boldsymbol{L} \boldsymbol{Y}=\boldsymbol{b}, \quad \boldsymbol{L}^{\mathrm{T}} \boldsymbol{X}=\boldsymbol{Y} LY=b,LTX=Y

{ y 1 = b 1 y i = b i − ∑ k = 1 i − 1 l i k y k , 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=bik=1i1likyk,i=2,3,,n

{ x n = y n / u n n x i = ( y i − ∑ k = i + 1 n u i k x k ) / u i i , 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=(yik=i+1nuikxk)/uii,i=n1,n2,,1.

由于上述线性方程组有大量开方运算,故改进Cholesky分解法。

改进的Cholesky分解法

A A A分解为:
A = L D L T A = LDL^T A=LDLT

Step 1. 计算 D 值

d j = a j j − ∑ k = 1 j − 1 l j k v j k , v j k = l j k d k , 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=ajjk=1j1ljkvjk,vjk=ljkdk,
在这里插入图片描述

Step 2. 计算第j列的L值

注意:L为单位下三角矩阵。即:L矩阵对角元素为1.
l i j = ( a i j − ∑ k = 1 j − 1 l i k v j k ) / d j , 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=(aijk=1j1likvjk)/dj,i=j+1,j+2,,n.
在这里插入图片描述

在完成分解后, 我们可分别求解下列系数矩。

L Y = b , D L T X = Y \boldsymbol{L Y}=\boldsymbol{b}, \quad \boldsymbol{D} \boldsymbol{L}^{\mathrm{T}} \boldsymbol{X}=\boldsymbol{Y} LY=b,DLTX=Y

追赶法

条件

  • 三对角非线性方程组
  • 对角占优( ∣ b i ∣ > = ∣ a i ∣ + ∣ c i ∣ |b_i|>=|a_i|+|c_i| bi>=ai+ci

( b 1 c 1 a 2 b 2 c 2 ⋱ ⋱ ⋱ a i b i c i ⋱ ⋱ ⋱ a n − 1 b n − 1 c n − 1 a n b n ) ( x 1 x 2 ⋮ x i ⋮ x n − 1 x n ) = ( d 1 d 2 ⋮ d i ⋮ d n − 1 d n ) \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) b1a2c1b2c2aibician1bn1ancn1bnx1x2xixn1xn=d1d2didn1dn

求解过程

  • 通过 Gauss 消元法将其化为系数矩阵为单位上三角形矩阵的方程组
    ( 1 β 1 1 β 2 ⋱ ⋱ 1 β n − 1 1 ) ( x 1 x 2 ⋮ x n − 1 x n ) = ( γ 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β21βn11x1x2xn1xn=γ1γ2γn1γn
  • 回代求得方程组的解。
    x n = γ n , x i = γ i − β i x i + 1 , i = n − 1 , n − 2 , ⋯   , 2 x_{n}=\gamma_{n}, \quad x_{i}=\gamma_{i}-\beta_{i} x_{i+1}, \quad i=n-1, n-2, \cdots, 2 xn=γn,xi=γiβixi+1,i=n1,n2,,2

误差分析

c o n d ( A ) = ∣ ∣ A − 1 ∣ ∣ ⋅ ∣ ∣ A ∣ ∣ cond(A) = ||A^{-1}||\cdot||A|| cond(A)=A1A刻画了方程组 A x = b Ax=b Ax=b的病态程度,及解对 A , b A,b A,b扰动的敏感程度。

性质

对可逆矩阵A,有:

  • c o n d ( A ) = c o n d ( A − 1 ) cond(A) = cond(A^{-1}) cond(A)=cond(A1)
  • c o n d ( k A ) = c o n d ( A ) cond(kA) = cond(A) cond(kA)=cond(A)
  • c o n d 2 ( U ) = 1 cond_2(U) = 1 cond2(U)=1
  • 等等
本文内容由网友自发贡献,版权归原作者所有,本站不承担相应法律责任。如您发现有涉嫌抄袭侵权的内容,请联系:hwhale#tublm.com(使用前将#替换为@)

计算方法--解线性方程组的直接法 的相关文章

  • 虚拟机网络适配器的三种模式详解及其配置

    VMWare中网络适配器的三种模式详解 关于虚拟机下Linux下ping www baidu com 出现 ping unknown host www baidu com问题的解决 有可能是因为网络适配器未正常配置 本文参照文章 https
  • Linux应用编程(系统信息与系统资源)

    在应用程序当中 有时往往需要去获取到一些系统相关的信息 譬如时间 日期 以及其它一些系统相关信息 本章将向大家介绍如何通过 Linux 系统调用或 C 库函数获取系统信息 譬如获取系统时间 日期以及设置系统时间 日期等 除此之外 还会向大家

随机推荐

  • Nacos介绍以及使用

    目录 一 概述 1 1 Nacos是什么 能干嘛 1 2 去哪下载 1 3 各个注册中心比较 二 Nacos作为服务注册中心 2 1 基于Nacos的服务提供者 2 2 基于Nacos的服务消费者 三 Nacos作为服务配置中心 3 1 N
  • 基于深度学习的海上雷达数据质量管控自动化技术

    文章作者Rune Gangeskar任职于Miros公司 目标是设计一套Wavex传感器系统 如何精准测量波浪 洋流 以及对水航速 并使用深度学习网络来自动辨识测量下取得的雷达数据 进一步提升Wavex系统的表现与可靠度 一 总体简介 对海
  • Java实现不规则软件版本号比较大小

    背景 最近由于需要比较两个版本号 从网上寻找的例子出现了问题 因此单独写一个不规则的版本号比较方法 代码 如果version1大于等于version2就返回true 可以根据自己需要进行调整 public static boolean co
  • 使用 Simulink 进行 STM32 编程

    目录 介绍 所需材料 步骤 1 在MATLAB中设置STM32 MAT软件路径 步骤 2 在STM32CubeMX中创建一个项目 步骤 3 配置时钟和 GPIO 引脚 步骤 4 项目经理并生成代码 步骤 5 在 Simulink 中创建模型
  • Java多线程——线程同步

    1 不安全的买票 不安全的买票 线程不安全 有负数或者多人买到同一张票 public class UnsafeBuyTicket public static void main String args BuyTicket buyTicket
  • $ 与 $$(*) 是什么?

    前言 今天在阅读一篇文章 小伙伴遇到这个问题说不想干前端了 一次Chrome翻译造成的玄学bug 主要内容是说翻译插件引发页面的react报错 文章中用 去获取报错页面所有的dom标签和正常无报错页面的dom标签 再进行对比来定位问题 具体
  • webform jquery ajax,ajax - Asp.net webform, jquery - Stack Overflow

    I m working on a project of estimation for my office I m using Ajax Autocomplete in Visual Studio 2017 The following cod
  • 第九章:构造器与垃圾收集器-对象的前世今生

    该系列文章系个人读书笔记及总结性内容 任何组织和个人不得转载进行商业活动 第九章 构造器与垃圾收集器 对象的前世今生 对象有生有死 必须对对象的生命周期负责 何时创建 何时销毁 不是消灭 只是放弃 由垃圾收集器 GC 将它蒸发掉 回收所占空
  • git命令总结

    1 git init 在当前目录下创建新的git仓库 2 git add filename 文件版本控制之前需要对这些文件进行追踪 对filename进行追踪 将文件添加进入缓存 3 git commit 提交更新 git commit a
  • IPv4 地址已耗尽,IPv6 涅槃重生:腾讯云IPv6改造综述

    引言 近日 全球 IPv4 地址正式耗尽的消息刷遍各大技术媒体 IPv6 再一次被推到人们面前 IP 作为网络世界的通行证 其重要性不言而喻 IPv4 地址枯竭 IPv6 作为IPv4地址枯竭的解决方案 其在中国的发展历程是怎样的 产品环环
  • 微信小程序开源项目精选

    本期为大家精选了码云上优秀的微信小程序开源项目 包括电商 博客 框架 建站系统 日常工具 图像识别等 希望能够给大家带来一点帮助 1 项目名称 微信电商小程序 作者 三三网络科技 项目简介 此项目是一套完整的电商系统 并且兼容各种电商场景可
  • 粒子滤波原理及其matlab仿真

    Github个人博客 https joeyos github io 粒子滤波原理及其matlab仿真 系统建模 粒子滤波算法不受线性高斯模型的约束 与卡尔曼滤波器一样 粒子滤波算法同样需要知道系统的模型 如果不知道系统的模型 也要想办法构建
  • Ubuntu 22.04 解决使用 .AppImage 文件方法

    AppImages 是一个文件系统 需要 FUSE 版本为 2 才能运行 但是 Ubuntu 22 04 的发行版本没有对其进行原始的配置的安装 重新安装并且配置即可 终端问题 未加载到 libfuse so 2 软连接 所以无法执行 Ap
  • C++对csv文件的读写

    include
  • 如何在线免费将MP4转换成MP3格式音乐

    MP4已经成为互联网上最流行的视频格式 我们从各种视频资源网站上下载到的视频文件大部分都是以MP4格式存储的 尤其是一些高品质的歌曲MV 为了达到在高压缩的前提下得到最好的质量 几乎都是mp4文件 但是如果你想直接在手机或者车上听这些歌曲
  • 【Mac】【Git】 全局配置 忽略 .DS_Store

    DS Store 是什么 Mac OS X 使用 DS Store 文件来存储文件夹特定的元数据信息 它们是在 Mac OS X Finder 访问的每个文件夹中创建的 甚至是网络宗卷和外部设备 文件夹级别的自定义存储在 DS Store
  • Vue锚链接(两种方法) scrollIntoView

    第一种 常见 锚链接 id和 href 结合起来 div style height 300px 第一 div div style height 300px 第二 div a href one 回到第一 a a href two 回到第二 a
  • 直流电源_滤波电路

    目录 前言 电容滤波电路 1 滤波原理 2 输出电压平均值 3 脉动系数 4 整流二极管的导通角 5 电容滤波电路的输出特性和滤波特性 倍压整流电路 其他形式的滤波电路 1 电感滤波电路 2 复式滤波电路 3 各种滤波电路的比较 前言 整流
  • 较好用的html5编译器

    一 BlueGriffon 基于 Firefox 渲染引擎的下一代 Web 和 EPUB 编辑器 链接 BlueGriffon 二 Aloha Editor 基于 所见即所得 的原则 HTML5 编辑器可以直接在门户网站上编辑网站 快速 简
  • 计算方法--解线性方程组的直接法

    文章目录 一 Gauss 消元法 1 顺序高斯消元法 总计算量 2 主元素高斯消元法 列主元素高斯消元法 3 高斯 约当 Gauss Jordan 消去法 总计算量 二 矩阵三角分解法 1 直接三角分解法 LU分解 Doolittle分解