• 解线性方程组python实现直接分解法(Doolittle,克劳特,追赶法)


    1. Doolittle分解

            Doolittle分解法是一种将矩阵分解为下三角矩阵L和上三角矩阵U的方法。它是LU分解的一种特殊形式,常用于线性方程组的求解和矩阵求逆等计算中。

            Doolittle分解的基本概念如下:

    1. LU分解:给定一个n×n的矩阵A,LU分解将其分解为一个下三角矩阵L和一个上三角矩阵U的乘积,即A = LU。

    2. 下三角矩阵L:下三角矩阵L的对角线元素都为1,其余元素满足L的上三角部分全为0。L的非零元素表示A矩阵中每个位置元素的消去倍数。

    3. 上三角矩阵U:上三角矩阵U的对角线元素和上三角部分非零,其余元素全为0。U的非零元素表示经过消元操作后得到的A矩阵的值。

    4. Doolittle分解:Doolittle分解是LU分解的一种形式,其中下三角矩阵L的对角线元素为1。通过Doolittle分解,可以将线性方程组的求解问题转化为先求解Ly = b,再求解Ux = y的过程。

            Doolittle分解的优点是不需要进行行交换操作,因此在某些特殊的情况下更加高效。它常被用于求解线性方程组、计算矩阵的逆以及求解矩阵特征值等应用中。

            请注意,Doolittle分解要求原始矩阵A的主子矩阵都非奇异,即存在唯一的分解。当遇到奇异矩阵或无法进行分解的情况时,需要使用其他方法进行处理。

    1. """
    2. @Time : 2023/10/19 0019 17:11
    3. @Auth : yeqc
    4. Doolittle分解
    5. """
    6. import numpy as np
    7. def doolittle_decomposition(A):
    8. n = len(A)
    9. L = np.zeros((n, n))
    10. U = np.zeros((n, n))
    11. for i in range(n):
    12. # 计算U的第一行
    13. U[0][i] = A[0][i]
    14. for i in range(n):
    15. # 计算L的第一列
    16. L[i][0] = A[i][0] / U[0][0]
    17. for i in range(1, n):
    18. for j in range(i, n):
    19. sum1 = 0.0
    20. for k in range(i):
    21. sum1 += L[i][k] * U[k][j]
    22. U[i][j] = A[i][j] - sum1
    23. for j in range(i + 1, n):
    24. sum2 = 0.0
    25. for k in range(i):
    26. sum2 += L[j][k] * U[k][i]
    27. L[j][i] = (A[j][i] - sum2) / U[i][i]
    28. return L, U
    29. def solve_linear_equations(A, b):
    30. L, U = doolittle_decomposition(A)
    31. n = len(A)
    32. y = np.zeros(n)
    33. x = np.zeros(n)
    34. # 解 Ly = b
    35. for i in range(n):
    36. sum1 = 0.0
    37. for j in range(i):
    38. sum1 += L[i][j] * y[j]
    39. y[i] = b[i] - sum1
    40. # 解 Ux = y
    41. for i in range(n-1, -1, -1):
    42. sum2 = 0.0
    43. for j in range(i+1, n):
    44. sum2 += U[i][j] * x[j]
    45. x[i] = (y[i] - sum2) / U[i][i]
    46. return x
    47. # 示例
    48. A = np.array([[9, -1, -1], [-1, 8, 0], [-1, 0, 9]])
    49. b = np.array([7, 7, 8])
    50. x = solve_linear_equations(A, b)
    51. print("解:", x)

    2. 克劳特分解

            克劳特分解法(Crout's method)是一种将矩阵分解为下三角矩阵L和上三角矩阵U的方法。它与Doolittle分解类似,但在求解过程中有一些区别。

            克劳特分解的基本概念如下:

    1. LU分解:给定一个n×n的矩阵A,LU分解将其分解为一个下三角矩阵L和一个上三角矩阵U的乘积,即A = LU。

    2. 下三角矩阵L:下三角矩阵L的对角线元素都为1,其余元素满足L的上三角部分全为0。L的非零元素表示A矩阵中每个位置元素的消去倍数。

    3. 上三角矩阵U:上三角矩阵U的对角线元素和上三角部分非零,其余元素全为0。U的非零元素表示经过消元操作后得到的A矩阵的值。

    4. 克劳特分解:克劳特分解是LU分解的一种形式,其中上三角矩阵U的对角线元素为1。通过克劳特分解,可以将线性方程组的求解问题转化为先求解Ly = b,再求解Ux = y的过程。

            克劳特分解与Doolittle分解的区别在于它将单位矩阵的元素放在上三角矩阵U的对角线上,而不是下三角矩阵L。这个特点使得在进行分解时,不需要对U的对角线元素进行除法运算,从而简化了计算过程。

            克劳特分解常用于求解线性方程组、计算矩阵的逆以及求解矩阵特征值等应用中。需要注意的是,克劳特分解要求原始矩阵A的主子矩阵都非奇异,即存在唯一的分解。当遇到奇异矩阵或无法进行分解的情况时,需要使用其他方法进行处理。        

    1. """
    2. @Time : 2023/10/19 0019 17:11
    3. @Auth : yeqc
    4. 克劳特分解
    5. """
    6. import numpy as np
    7. def crout_decomposition(A):
    8. n = len(A)
    9. L = np.zeros((n, n))
    10. U = np.zeros((n, n))
    11. for i in range(n):
    12. # 计算U的第一行
    13. U[0][i] = A[0][i]
    14. for i in range(n):
    15. # 计算L的第一列(除了对角线元素为1)
    16. L[i][0] = A[i][0] / U[0][0]
    17. for i in range(1, n):
    18. for j in range(i, n):
    19. sum1 = 0.0
    20. for k in range(i):
    21. sum1 += L[j][k] * U[k][i]
    22. U[j][i] = A[j][i] - sum1
    23. for j in range(i, n):
    24. sum2 = 0.0
    25. for k in range(i):
    26. sum2 += L[i][k] * U[k][j]
    27. L[i][j] = (A[i][j] - sum2) / U[i][i]
    28. return L, U
    29. def solve_linear_equations(A, b):
    30. L, U = crout_decomposition(A)
    31. n = len(A)
    32. y = np.zeros(n)
    33. x = np.zeros(n)
    34. # 解 Ly = b
    35. for i in range(n):
    36. sum1 = 0.0
    37. for j in range(i):
    38. sum1 += L[i][j] * y[j]
    39. y[i] = b[i] - sum1 / L[i][i]
    40. # 解 Ux = y
    41. for i in range(n - 1, -1, -1):
    42. sum2 = 0.0
    43. for j in range(i + 1, n):
    44. sum2 += U[i][j] * x[j]
    45. x[i] = (y[i] - sum2) / U[i][i]
    46. return x
    47. # 示例
    48. A = np.array([[9, -1, -1], [-1, 8, 0], [-1, 0, 9]])
    49. b = np.array([7, 7, 8])
    50. x = solve_linear_equations(A, b)
    51. print("解:", x)

    3. 三对角阵的追赶法

            三对角矩阵是指只有主对角线及其两侧的一条对角线上有非零元素,其余位置都是0的方阵。

            三对角矩阵在数值计算中非常重要,因为许多微分方程离散化后可以转化为三对角矩阵的求解问题。使用标准的高斯消元法解三对角矩阵是比较困难的,因此需要一种特殊的算法来高效地求解。

            追赶法(也称为Thomas算法)是一种用于解三对角线性方程组的直接求解方法。它的主要思想是通过高斯消元的思想将系数矩阵化为上三角矩阵,再利用回代的方法求解方程组。由于三对角矩阵具有特殊的结构,追赶法只需要进行简单的数学运算即可高效地求解线性方程组,时间复杂度为O(n)。

    1. """
    2. @Time : 2023/10/19 0019 17:11
    3. @Auth : yeqc
    4. 追赶法
    5. """
    6. import numpy as np
    7. def solve_tridiagonal_system(a, b, c, d):
    8. n = len(d)
    9. x = np.zeros(n)
    10. # 将系数矩阵化为上三角矩阵
    11. for i in range(1, n):
    12. m = a[i] / b[i - 1]
    13. b[i] -= m * c[i - 1]
    14. d[i] -= m * d[i - 1]
    15. # 回代求解
    16. x[n - 1] = d[n - 1] / b[n - 1]
    17. for i in range(n - 2, -1, -1):
    18. x[i] = (d[i] - c[i] * x[i + 1]) / b[i]
    19. return x
    20. # 示例使用
    21. a = np.array([0, 2, 2]) # 下对角线元素
    22. b = np.array([4, 4, 4]) # 主对角线元素
    23. c = np.array([2, 2, 0]) # 上对角线元素
    24. d = np.array([6, 12, 10]) # 右侧常数向量
    25. x = solve_tridiagonal_system(a, b, c, d)
    26. print("解为:", x)

  • 相关阅读:
    【建议收藏】Kafka 面试连环炮, 看看你能撑到哪一步?
    c++中的模板(8) -- 模板和友元函数
    技术学习:Python |欲先善其事,必先利其器(JSON)二
    (分布式缓存)Redis分片集群
    二百三十三、Flume——Flume采集JSON文件到Kafka,再用Flume采集Kafka数据到HDFS中
    C++中istream_iterator和ostream_iterator的源码分析
    隔离出来的“陋室铭”
    学习笔记:机器学习之支持向量机(七、求核函数)
    小白备战大厂算法笔试(七)——图
    Python批量清除目录结构保留文件
  • 原文地址:https://blog.csdn.net/anxinbuxinan/article/details/133957830