列主元高斯消去法数学原理及超级纯手工Python实现

一、引言

从高斯消去法,我们看到还是有缺陷,高斯消去法中,当对角元素=0时,消去无法进行,当对角元素很小的时候,导致其他它元素数量级严重增长和舍入误差扩散1,使结果不可靠。因此引出了列主元高斯消去法。

解决办法:交换行,避免小主元。

每次消元前,在要消去未知数的系数中找到绝对值最大的系数作为主元,通过方程对换将其换到对角线上,然后进行消元。

二、python代码实现

第一步:找主元交换行:

n=A.shape[0]
for i in range(n-1):
        max_=0
        for j in range(i,n):
            if np.abs(A[j,i])>max_:
                max_=np.abs(A[j,i])
                maxn=j
        ##将第i行与第maxn行互换
        exchange=A[i].copy()
        A[i]=A[maxn]
        A[maxn]=exchange

结果如图所示:

原矩阵:

换完后:

 这块一定要注意!!!!!!一定得copy()不然交换最后一行不会改变。还有一点在找最大主元的时候不建议掉函数,多少还是有问题的!可以思考下面这个图片,但是在消元中出现的。

剩下的消元和回代和消去法思想一样,但是我们要注意A是增广矩阵了!!!

(二)函数代码:

def sloveeq(A,a):#A是增广矩阵
    if a.shape[0]!=a.shape[1]:
        print('系数矩阵不是方阵')
        sys.exit(0)
    n=A.shape[0]
    for i in range(n-1):
        max_=0
        for j in range(i,n):
            if np.abs(A[j,i])>max_:
                max_=np.abs(A[j,i])
                maxn=j
        ##将第i行与第maxn行互换
        exchange=A[i].copy()
        A[i]=A[maxn]
        A[maxn]=exchange
        #消元
    for k in range(i+1,n):
            q=A[k,i]/A[i,i]
            A[k,:]= A[k,:]-q*A[i,:]
    #回代
    b=A[:,n]
    b[n-1]=A[n-1,n]/A[n-1,n-1]
    for i in range(n-2,-1,-1):
        b[i]=(b[i]-np.dot(A[i,i+1:n],b[i+1:n]))/A[i,i]
    return b

三、实验:

a=np.array([[0.001,2,3],[-1,3.712,4.623],[-2,1.072,5.643]])
c=np.array([1.0,2.0,3.0])
A= np.insert(a, 3, values=c, axis=1)
X=sloveeq(A,a)
X

结果如下:

(-0.49039646, -0.05103518,  0.3675202)
  • 3
    点赞
  • 26
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
Python中的列主元高斯消去法是一种用于解线性方程组的数值方法。它通过将线性方程组表示为增广矩阵的形式,并利用矩阵的行变换来消去未知数,最终得到方程组的解。 以下是Python实现列主元高斯消去法的基本步骤: 1. 构建增广矩阵:将线性方程组的系数矩阵和常数向量合并成一个增广矩阵。 2. 选取列主元:在每一次消元操作中,选择当前中绝对值最大的元素作为主元。 3. 行交换:将含有主元的行与当前操作行交换,确保主元位于当前操作行的第一个位置。 4. 消元计算:通过对当前操作行进行线性组合,将当前下方的所有元素消为零。 5. 重复步骤2-4,直到得到上三角矩阵形式。 6. 回代求解:从最后一行开始,通过回代计算出未知数的值。 下面是一个简单的Python代码示例,演示了如何使用列主元高斯消去法解线性方程组: ```python import numpy as np def gaussian_elimination(A, b): n = len(b) Ab = np.concatenate((A, b.reshape(n, 1)), axis=1) for i in range(n): # 选取列主元 max_index = np.argmax(np.abs(Ab[i:, i])) + i Ab[[i, max_index]] = Ab[[max_index, i]] # 消元计算 for j in range(i+1, n): ratio = Ab[j, i] / Ab[i, i] Ab[j] -= ratio * Ab[i] # 回代求解 x = np.zeros(n) for i in range(n-1, -1, -1): x[i] = (Ab[i, -1] - np.dot(Ab[i, :-1], x)) / Ab[i, i] return x # 示例方程组的系数矩阵和常数向量 A = np.array([[2, 1, -1], [-3, -1, 2], [-2, 1, 2]]) b = np.array([8, -11, -3]) # 调用函数求解方程组 solution = gaussian_elimination(A, b) print("方程组的解为:", solution) ``` 这段代码使用了NumPy库来进行矩阵运算,首先定义了一个`gaussian_elimination`函数,接受系数矩阵A和常数向量b作为输入,返回方程组的解。在函数内部,通过循环实现列主元高斯消去法的各个步骤,最后得到方程组的解并打印输出。

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值