我想用python的numpy和sympy库实现以下功能，首先已知一个方阵A，其中的元素都是具体的数值，接着有一个列矩阵b，它的行数和A相等，同时其中的元素是符号比如（a1，a2……）这样子的，我想在把这两个矩阵横向连接，形成一个增广矩阵。接着我想要通过初等变换的形式将这个矩阵中原本是方阵的部分转化为一个单位矩阵，当然每一次的初等变换中列矩阵b也要跟着变化。最后我想要得到那个新的带有符号的列矩阵b。请注意由于列矩阵b中原本的元素都是符号，因此在做初等变换的时候就用sympy库来对其进行变换在得到new_b这个列向量之后，我想要再把其中每个元素对于原本的b向量中的每个元素的系数提取出来，形成一个新的方阵，请问代码要怎么写

可以按照以下步骤来实现：

1. 导入必要的库

import numpy as np
import sympy as sp

2. 创建方阵A和列向量b

A = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
b = sp.Matrix(sp.symbols('a1:4'))

3. 将A和b横向连接形成增广矩阵

augmented_matrix = np.hstack((A, np.array(b).reshape((-1, 1))))

4. 对增广矩阵进行初等变换，将A部分转化为单位矩阵

n_rows, n_cols = augmented_matrix.shape
for i in range(n_rows):
    if i >= n_cols - 1:
        break
    pivot = augmented_matrix[i, i]
    if pivot == 0:
        continue
    augmented_matrix[i, :] /= pivot
    for j in range(i+1, n_rows):
        factor = augmented_matrix[j, i] / augmented_matrix[i, i]
        augmented_matrix[j, :] -= factor * augmented_matrix[i, :]

5. 对列矩阵b进行相同的初等变换

for i in range(n_rows):
    if i >= n_cols - 1:
        break
    pivot = augmented_matrix[i, i]
    if pivot == 0:
        continue
    b[i] /= pivot
    for j in range(i+1, n_rows):
        factor = augmented_matrix[j, i] / augmented_matrix[i, i]
        b[j] -= factor * b[i]

6. 从new_b中提取系数形成新的方阵

new_b = sp.Matrix(b)
coeff_matrix = sp.Matrix(np.zeros((n_rows, n_rows)))
for i in range(n_rows):
    for j in range(n_cols-1, n_rows):
        coeff_matrix[i, j-n_cols+1] = new_b[i].coeff(sp.Symbol('a{}'.format(j-n_cols+2)))

完整代码如下：

import numpy as np
import sympy as sp

A = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
b = sp.Matrix(sp.symbols('a1:4'))

augmented_matrix = np.hstack((A, np.array(b).reshape((-1, 1))))

n_rows, n_cols = augmented_matrix.shape
for i in range(n_rows):
    if i >= n_cols - 1:
        break
    pivot = augmented_matrix[i, i]
    if pivot == 0:
        continue
    augmented_matrix[i, :] /= pivot
    for j in range(i+1, n_rows):
        factor = augmented_matrix[j, i] / augmented_matrix[i, i]
        augmented_matrix[j, :] -= factor * augmented_matrix[i, :]

for i in range(n_rows):
    if i >= n_cols - 1:
        break
    pivot = augmented_matrix[i, i]
    if pivot == 0:
        continue
    b[i] /= pivot
    for j in range(i+1, n_rows):
        factor = augmented_matrix[j, i] / augmented_matrix[i, i]
        b[j] -= factor * b[i]

new_b = sp.Matrix(b)
coeff_matrix = sp.Matrix(np.zeros((n_rows, n_rows)))
for i in range(n_rows):
    for j in range(n_cols-1, n_rows):
        coeff_matrix[i, j-n_cols+1] = new_b[i].coeff(sp.Symbol('a{}'.format(j-n_cols+2)))

print(coeff_matrix)