高斯消元法求解方程组(要有python基础和线性代数的基础)
本人这学期开的是一门数值分析课,老师要求用python写出实现某些算法的代码,当遇到第一个高斯消元法,印象深刻的是这个编程与以往的编程不一样,从这几点来讲:首先,我是在上这门课之前就将python的基础语法过了一遍,但是第一次看着老师的代码,发现下标标索引如此之多的,心里乱哄哄的;其次,有着c语言的编程基础,我一次做多看见了三个for循环,再加上低一点,我整个人就懵了;然后再看到网上没有详细的过程,这样我也很难理解。
接下来我就想将整个方法的动画演示,数学公式,有了这些再接下来对比着代码看就会比较轻松。
目录
一.基本理论
(1)基本原理
二.思路讲解
1)消元
(1)基本思路
(2)举个栗子
(3)规律总结:
(4)代码理解
2)回代过程
(1)基本思路:
(2)代码实现
三.总体代码实现
四.总结
1.消元的思路:我并没有采用一行一行的消元,而是采用一列一列的消元,一行一行的我也试过,感觉这样比较繁琐。我会在接下来继续试一试如何去写出一行一行消元的代码。
2.可能对于python知识的欠缺,我这个搞了三天,代码还是听我们老师讲了一遍(只讲了一半),真的很痛苦,直到昨天我理解了代码的思路,然后也想着用两种方法(另外一种就是按行进行消元,但是却是以失败告终),这也是我迟迟没有更新的原因,谢谢大家!这真是一个联系编程能力的好课程。
一.基本理论
(1)基本原理
高斯消元法主要是针对有解的并且解是唯一的矩阵(当det(A)!=0时),通过先构造一个增广矩阵,然后对其进行操作行初等变化,得到上三角矩阵(主对角线元素以下的所有元素都为0,主对角线上的元素可以为1,也可以不为1)。然后从下往上依次求解出所有的解,即第一步求解出x1,第二步运用第一步求解的结果求出第二个解x2,第三步运用前面求出的解再求出第三个解,……,第n步运用前面求出的n-1个解求出第n个解。主要比较难理解的还是消元过程。
二.思路讲解
矩阵的变化主要过程有:
1)消元
(1)基本思路
消去主对角线下面的所有元素,按照从上到下、从左到右的顺序。
增广矩阵 - > 上三角矩阵:下面依次就是增广矩阵和上三角矩阵
那么这一部分之间是如何用python实现的,也就说需要将表达式写出来,通过写一个简单的例子,我们先来看一下:
(2)举个栗子
我们发现上面的增广矩阵中对角线一下2,4,-9根据上三角矩阵的概念都是要消去的,而我们的
思想是先消去2(利用第一行),再消去4(利用第一行),再消去-9(利用第二行)。我们会发现
这个是有一定规律的。我们先看一下这个增广矩阵通过初等变化得到上三角矩阵的计算过程。
1.消去2,对于人来说非常简单,但是对于计算机的话我们需要从计算机的角度考虑问题:
我们想想当第一行乘以多少的时候加在第二行会将第二行的要第一个元素消去,现在的问题就转化
成求乘的那个数?设这个数为x(方便理解,当然你也可以不用看)即满足1*x+2=0,也就说
x=-2/1,写成索引的形式也就是x=-a[1][0] / a[0][0]。先理解这一步是非常重要的。然后就会得到新的
第二行用python代码写就是(a[0]为第一行,a[0][0]为第一行第一列的第一个元素,这些是提前必须
知道的):a[1]=a[0]*x+a[1] 也就是 a[1]=a[0]* -a[1][0] / a[0][0] +a[1] 左边的新的第二行,右边
的是原来的第二行。
2. 消去4,4处于第三行第一列,因为上一步骤已经将第二行的第一列元素消去了,所以只能用第一
行第一列的数据消去4,和上面一样,设x找4与第一行第一列元素的代数关系,和第一步一样
a[2]= a[0]*x+a[2] 即,a[2]= a[0]* -a[2][0] / a[0][0] +a[2]
3.消去-9,利用上面求出的新行,即用第二行。因为此时如果用第一行或者其他行,这样会导致上
一次消去4的位置会重新出现一个新的数。
(3)规律总结:
这一步会方便我们在接下来的工作中推导出python代码:
从上面我们可以发现,根据numpy的特点和python代码的特点,我们按列进行消元会比较方便:就是依次先将主对角线上第一个元素下面的数值(有n-1个)消为0,然后是消主对角线上第二个元素下面的元素(有n-2个),……,然后消到倒数第二个(因为倒数第一个的下面没有元素可以消了)
我们会发现在这里会遇到两个for循环,第一个for循环 for col in range(len(m[0])):
for 循环 for row in range(col+1,len(m)):一个循环代表消去哪一行的零。必须至少有两个循环,一个控制列,一个控制行,在这里停顿一下,思考一下:
下面看看我写的代码,在代码里面有详细的介绍。
for col in range(len(m)):#按列进行消元,一次将一个矩阵主对角线上的一个主角元素下边的所有元素消为0for row in range(col+1,len(m)):#越往下循环的次数越来越少,第一个主角元素下面共有n-1,也就是说要消去n-1个元素,第二个主角元素需要消去n-2个元素,最后一个主角元素下面没有元素,所以就不需要消去一些元素,即刚好遍历到n-1r = [(rowValue*(-(m[row][col]/m[col][col]))) for rowValue in m[col]] #rowValue表示某一行的元素,对于内for循环来说,rowvalue表示第一行的数据,用这一行要将下面的所有行的第一个元素消去,当然下面每一行的第一个数前面的系数不可能是一样的,因此,-(m[row][col]/m[col][col])表示每一行与第一行的系数关系,每一步都需要计算,然后计算出可以消去的行,再在下面的列表表达式中计算。m[row] = [sum(pair) for pair in zip(m[row], r)]#得到新的row行 m[r]:是经过没有消零以前的行与r加法的结果,这是列表表达式的基础用法。
(4)代码理解
这是两个for循环的嵌套,这是我经过调试发现的:第一个for遍历到0时,m[0]表示一整行的数据,然后到第二个循环时,要消去的有n-1个(主对角线第一个元素下面有n-1个数)。看看我上面举得例子方便理解哈,再对比代码去看看。
2)回代过程
(1)基本思路:
进行完消元过程之后我们得到了一个上三角矩阵,就是下面这种上面元素最多,下面只有一个元素了,根据我们回代的过程,我们是从下面回代的,从计算机的角度分析看:会发现,将矩阵reverse()方法一下计算起来会比较方便,计算出结果之后在reverse()方法一下,等于说计算出来的值顺序不会发生改变。
(2)代码实现
ans=[]#初始化一个矩阵m = list(m)m.reverse()for sol in range(len(m)):if sol == 0:ans.append(m[sol][-1]/m[sol][-2])#倒数第一个数除以倒数第二个数,这是第一行的情况#ans里面放方程的根else:inner = 0#这是第一行下面的情况#for x in range(sol):#inner += (ans[x]*m[sol][-2-x])#将求出的所有解依次乘以其下面的元素值# the equation is now reduced to ax+b=c form# solve with (c-b)/aans.append((m[sol][-1]-inner)/m[sol][-sol-2])#右边的值减去上一个的结果然后除以x的系数(这一行的倒数第三个元素)
对于求出的第一个解时不用考虑前面的结果;但是当计算第二个解时,有解线性方程组的知识可以找到,需要考虑第一个解(看我上面的代码实现过程);当计算第三个解时,需要考虑前面两个解的情况。
三.总体代码实现
"""
作者:小翟同学
日期:2022年09月19日
"""
import numpy as np
m = np.array([[1.0,-1.0,3.0,1.0], [2.0, -4.0, 6.0,4.0], [4.0, -9.0, 2.0,1.0]])def myGauss(m):# eliminate columnsfor col in range(len(m[0])):#按列进行消元,一次将一个矩阵主对角线上的一个主角元素下边的所有元素消为0for row in range(col+1,len(m)):#越往下循环的次数越来越少,第一个主角元素下面共有n-1,也就是说要消去n-1个元素,第二个主角元素需要消去n-2个元素,最后一个主角元素下面没有元素,所以就不需要消去一些元素,即刚好遍历到n-1r = [(rowValue*(-(m[row][col]/m[col][col]))) for rowValue in m[col]] #rowValue表示每一行的元素,r的意思就是,m[col][col]( 吧版本信息表)除以要消去的的元素的系数,得到一个第一行与要消去的这一行的关系比例,然后乘在第一行# m[row] = [sum(pair) for pair in zip(m[row], r)]#得到新的row行为r(上面对第一行处理过的可以,可以得到一个改行下面的一个元素的额为0)+旧的m[row]#举个例子:m=[1,2,3]# r=[4,5,6]# m=[sum(i)for i in zip(m,r)]#运行结果 m# [5, 7, 9]# 在经过上面的过程之后 我们就可以得到一个主对角线元素为1,并且为上三角矩阵的行列式,#接下来就需要进行回代ans=[]#初始化一个矩阵m = list(m)m.reverse()for sol in range(len(m)):if sol == 0:ans.append(m[sol][-1]/m[sol][-2])#倒数第一个数除以倒数第二个数,这是第一行的情况#ans里面放方程的根else:inner = 0#这是第一行下面的情况#for x in range(sol):#inner += (ans[x]*m[sol][-2-x])#将求出的所有解依次乘以其下面的元素值# the equation is now reduced to ax+b=c form# solve with (c-b)/aans.append((m[sol][-1]-inner)/m[sol][-sol-2])#右边的值减去上一个的结果然后除以x的系数(这一行的倒数第三个元素)ans.reverse()return ansprint(myGauss(m))
四.总结
1.消元的思路:我并没有采用一行一行的消元,而是采用一列一列的消元,一行一行的我也试过,感觉这样比较繁琐。我会在接下来继续试一试如何去写出一行一行消元的代码。
2.可能对于python知识的欠缺,我这个搞了三天,代码还是听我们老师讲了一遍(只讲了一半),真的很痛苦,直到昨天我理解了代码的思路,然后也想着用两种方法(另外一种就是按行进行消元,但是却是以失败告终),这也是我迟迟没有更新的原因,谢谢大家!这真是一个联系编程能力的好课程。
高斯消元法求解方程组(要有python基础和线性代数的基础)相关推荐
- 【机器学习算法专题(蓄力计划)】五、机器学习中的线性代数的基础操作
文章目录 线性代数的基础操作 简单生成一个行向量 行向量的转置 (无效做法) 行向量的转置 直接生成一个列向量 (二维数组) 向量的内积 向量的外积 生成一个零矩阵 生成一个对角矩阵 生成一个单位矩阵 ...
- 高斯消元法求解方程组
引子 高斯消元法简介 引例 求解过程 编程实现高斯消元法C 1. 引子 1. 高斯消元法简介 数学上,高斯消元法(Gaussian Elimination),是线性代数规划中的一个算法,可用来为线性方 ...
- 高斯消元法求解线性方程组(附python代码)
输入:a是m×n的系数矩阵,b是m×1的(列)向量. 输出:方程组的通解. 用高斯消元法(行化简法)解线性方程组步骤 1.构造方程组的增广矩阵 2.从最左边列往右,使用行化简算法把增广矩阵化为阶梯形, ...
- 如何用python求解方程组_用Python的Numpy求解线性方程组
在本文中,您将看到如何使用Python的Numpy库解决线性方程组. 什么是线性方程组? 维基百科将线性方程组定义为:在数学中,线性方程组(或线性系统)是两个或多个涉及同一组变量的线性方程的集合. 解 ...
- 高斯消元法解矩阵方程c语言,c++高斯消元法求解线性方程组
//高斯消元法求解方程组 #include #include using namespace std; #define MaxNum 10 int array[MaxNum][MaxNum] = { ...
- python解多元多次方程组_Python求解多重或非线性方程,python,多元,多次,方程组,线性方程组...
背景: 如何使用python求解多元多次方程组或者非线性方程组. 原创内容,转载注明出处!请勿用于商业用途! (上篇用python拟合2019nCov感染人数的文章被不少博主转载了,发的比较早,不少博 ...
- python线性方程组求解_python求解方程组的三种方法
python求解方程组的三种方法: Numpy求解方程组x + 2y = 3 4x + 5y = 6 当然我们可以手动写出解析解,然后写一个函数来求解,这实际上只是用 Python 来单纯做" ...
- Python 克莱姆法则求解方程组
克莱姆法则求解方程组 请用克莱姆法则解下面的线性方程2x2系统: 编写程序,用户输入数字a.b.c.d.e和f,然后显示x和y的结果.如果ad-bc为零,呈现"The equation ha ...
- 科学计算:Python VS. MATLAB(3)----线性代数基础
科学计算:Python VS. MATLAB(3)----线性代数基础 按:在介绍工具之前先对理论基础进行必要的回顾是很必要的.没有理论的基础,讲再多的应用都是空中楼阁.本文主要设涉及线性代数和矩阵论 ...
最新文章
- ccnp bcmsn 高级交换学习笔记2
- C# json解析字符串总是多出双引号_在JavaScript应用中将CSV转换为JSON
- 微信小程序tabBar不显示的问题描述解决
- KVM中I/O设备直接分配和SR-IOV(十六)
- 『树上匹配 树形dp』
- 《ACM国际大学生程序设计竞赛题解Ⅰ》——模拟题
- Siamese Network理解
- oracle 取系统当前年份_Oracle 获取当前日期及日期格式
- 单片机c语言模块化实例程序设计,单片机C语言模块化设计
- linux ntp软件下载,Linux_Linux时区同步问题(安装ntp软件过程),下载了一个windows的NTP服务程序 - phpStudy...
- python 字符串 f_Python格式化字符串(f,F,format,%)
- Oracle数据库merge into的使用,存在则更新,不存在则插入
- oracle共享内存段手工清理
- Java EnumMap工作原理及实现
- Unity3D实现AB包加载资源
- Xshell上传文件方法
- DSP CCS12.00 芯片:TMS320F28335 外部中断 XINT1, 和映射区域的 k1 -- k4 按键的功能实现
- excel npoi 连接_C#操作Excel(NPOI)
- 走出软件作坊:三五个人十来条枪 如何成为开发正规军 链接[收藏]
- 上海计算机应用基础考试培训班,上海市计算机一级考试辅导
热门文章
- 2016第19本:至关重要的关系
- Day05java-继承
- php dir opendir,php中目录操作opendir()、readdir()及scandir()用法示例
- 迪文总线摄像头倒车影像应用方案
- 12Python爬虫---Fiddler抓包工具使用
- Vue项目目录结构介绍讲解
- 关于项目启动会和项目开工会议的区别
- 数学计算机代码,GeoGebra(数学图形计算器)(示例代码)
- win7装xp双系统_win7配置最低要求是什么
- android rar工具,安卓解压缩软件-ZArchiver Pro 解压工具下载v0.8.5 安卓版-安卓解压神器西西软件下载...