本人这学期开的是一门数值分析课，老师要求用python写出实现某些算法的代码，当遇到第一个高斯消元法，印象深刻的是这个编程与以往的编程不一样，从这几点来讲：首先，我是在上这门课之前就将python的基础语法过了一遍，但是第一次看着老师的代码，发现下标标索引如此之多的，心里乱哄哄的；其次，有着c语言的编程基础，我一次做多看见了三个for循环，再加上低一点，我整个人就懵了；然后再看到网上没有详细的过程，这样我也很难理解。
接下来我就想将整个方法的动画演示，数学公式，有了这些再接下来对比着代码看就会比较轻松。
目录
1.消元的思路：我并没有采用一行一行的消元，而是采用一列一列的消元，一行一行的我也试过，感觉这样比较繁琐。我会在接下来继续试一试如何去写出一行一行消元的代码。
一.基本理论
（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)):#按列进行消元，一次将一个矩阵主对角线上的一个主角元素下边的所有元素消为0
for row in range(col+1,len(m)):#越往下循环的次数越来越少，第一个主角元素下面共有n-1,也就是说要消去n-1个元素,第二个主角元素需要消去n-2个元素，最后一个主角元素下面没有元素，所以就不需要消去一些元素，即刚好遍历到n-1
r = [(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)/a
ans.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 columns
for col in range(len(m[0])):#按列进行消元，一次将一个矩阵主对角线上的一个主角元素下边的所有元素消为0
for row in range(col+1,len(m)):#越往下循环的次数越来越少，第一个主角元素下面共有n-1,也就是说要消去n-1个元素,第二个主角元素需要消去n-2个元素，最后一个主角元素下面没有元素，所以就不需要消去一些元素，即刚好遍历到n-1
r = [(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)/a
ans.append((m[sol][-1]-inner)/m[sol][-sol-2])#右边的值减去上一个的结果然后除以x的系数（这一行的倒数第三个元素）
ans.reverse()
return ans
print(myGauss(m))
1108
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
