线性代数笔记
第1章 线性代数中的线性方程组
1.1线性方程组
线性方程组的定义
解集，三种情况：无解、唯一解、无穷多解
线性方程组可以用矩阵表示：系数矩阵、增广矩阵
解线性方程组思路：等价变换（化简）
化简线性方程组的三种基本变换：
- 这里用method,i,j,k表示，method<0:行变换;method>0:列变换。
- |method|=1表示用k乘某一行(列);|method|=2表示两行(列)互换;
- |method|=3表示将i行(列)乘以数k后加到j行(列)上。i,j均从0开始。
线性方程组的变换对应于增广矩阵的初等变换：
- 倍加变换、对换变换、倍乘变换
2个线性方程组的增广矩阵是行等价的，则它们有相同的解集（想一下为什么）
线性方程组的2个基本问题：是否相容，有解的话解是否唯一
练习题
3.数组(3,4,-2)是否为下列方程组的解？
为解线性方程组，编写代码如下，可通过transform_matrix函数对增广矩阵进行初等变换。解题步骤在代码后面。
import numpy as np
def transform_matrix(matrix,method=0,i=0,j=0,k=1):
#method=0:不作变换;method<0:行变换;method>0:列变换
#|method|=1表示用k乘某一行(列);|method|=2表示两行(列)互换;|method|=3表示将i行(列)乘以数k后加到j行(列)上
#i,j均从0开始
#1.生成初等矩阵
if method==0:
return matrix
elif method<0:
if i<0 or i>=matrix.shape[0]:
print("error i")
if j<0 or j>=matrix.shape[0]:
print("error j")
elif method>0:
if i<0 or i>=matrix.shape[1]:
print("error i")
if j<0 or j>=matrix.shape[1]:
print("error j")
d=matrix.shape[0] if method<0 else matrix.shape[1]
A=np.eye(d)
if abs(method)==1:
A[i,i]=k
elif abs(method)==2:
A[i,j]=1
A[j,i]=1
A[i,i]=0
A[j,j]=0
elif method==-3:
A[j,i]=k
elif method==3:
A[i,j]=k
else:
print("method error")
return matrix
if method<0:
ret=np.matmul(A,matrix)
return ret
elif method>0:
ret=np.matmul(matrix,A)
return ret
else:
return matrix
a=np.array([[1,2,3,4],[5,6,7,8],[9,10,11,12]]);a
array([[ 1, 2, 3, 4],
[ 5, 6, 7, 8],
[ 9, 10, 11, 12]])
b=transform_matrix(a,-1,2,0,10.5);b
array([[ 1. , 2. , 3. , 4. ],
[ 5. , 6. , 7. , 8. ],
[ 94.5, 105. , 115.5, 126. ]])
c=transform_matrix(b,-2,1,2,0);c
array([[ 1. , 2. , 3. , 4. ],
[ 94.5, 105. , 115.5, 126. ],
[ 5. , 6. , 7. , 8. ]])
d=transform_matrix(c,-3,0,2,10);d
array([[ 1. , 2. , 3. , 4. ],
[ 94.5, 105. , 115.5, 126. ],
[ 15. , 26. , 37. , 48. ]])
练习题3解题步骤：
a=np.array([[5,-1,2,7],[-2,6,9,0],[-7,5,-3,-7]]);a
array([[ 5, -1, 2, 7],
[-2, 6, 9, 0],
[-7, 5, -3, -7]])
b=transform_matrix(a,-3,0,1,2/5);b
array([[ 5. , -1. , 2. , 7. ],
[ 0. , 5.6, 9.8, 2.8],
[-7. , 5. , -3. , -7. ]])
c=transform_matrix(b,-3,0,2,7/5);c
array([[ 5. , -1. , 2. , 7. ],
[ 0. , 5.6, 9.8, 2.8],
[ 0. , 3.6, -0.2, 2.8]])
d=transform_matrix(c,-3,1,2,-3.6/5.6);d
array([[ 5. , -1. , 2. , 7. ],
[ 0. , 5.6, 9.8, 2.8],
[ 0. , 0. , -6.5, 1. ]])
e=transform_matrix(d,-3,2,1,9.8/6.5);e
array([[ 5.00000000e+00, -1.00000000e+00, 2.00000000e+00,
7.00000000e+00],
[ 0.00000000e+00, 5.60000000e+00, -1.89421128e-15,
4.30769231e+00],
[ 0.00000000e+00, 0.00000000e+00, -6.50000000e+00,
1.00000000e+00]])
出问题了，是舍入误差导致的，想了半天没找到从程序上根本解决的办法，只有迂回解决了：
e=transform_matrix(d,-1,1,0,65);e
array([[ 5. , -1. , 2. , 7. ],
[ 0. , 364. , 637. , 182. ],
[ 0. , 0. , -6.5, 1. ]])
f=transform_matrix(e,-1,2,0,98);f
array([[ 5., -1., 2., 7.],
[ 0., 364., 637., 182.],
[ 0., 0., -637., 98.]])
g=transform_matrix(f,-3,2,1,1);g
array([[ 5.00000000e+00, -1.00000000e+00, 2.00000000e+00,
7.00000000e+00],
[ 0.00000000e+00, 3.64000000e+02, -1.13686838e-13,
2.80000000e+02],
[ 0.00000000e+00, 0.00000000e+00, -6.37000000e+02,
9.80000000e+01]])
还是不行，通过调试发现从d=transform_matrix(c,-3,1,2,-3.6/5.6)的执行过程中已产生了误差，d[2,2]=-6.500000000000001。百度关键词numpy 舍入误差，找到《numpy矩阵求逆舍入误差(numpy matrix inversion rounding errors)》，通过在程序中加入np.set_printoptions(suppress=True)即可改进显示效果。
np.set_printoptions(suppress=True)
g
array([[ 5., -1., 2., 7.],
[ 0., 364., -0., 280.],
[ 0., 0., -637., 98.]])
h=transform_matrix(g,-3,2,0,2/637);h
array([[ 5. , -1. , -0. , 7.30769231],
[ 0. , 364. , -0. , 280. ],
[ 0. , 0. , -637. , 98. ]])
i=transform_matrix(h,-3,1,0,1/364);i
array([[ 5. , 0. , -0. , 8.07692308],
[ 0. , 364. , -0. , 280. ],
[ 0. , 0. , -637. , 98. ]])
j=transform_matrix(i,-1,0,0,1/5);j
array([[ 1. , 0. , -0. , 1.61538462],
[ 0. , 364. , -0. , 280. ],
[ 0. , 0. , -637. , 98. ]])
k=transform_matrix(j,-1,1,0,1/364);k
array([[ 1. , 0. , -0. , 1.61538462],
[ 0. , 1. , -0. , 0.76923077],
[ 0. , 0. , -637. , 98. ]])
l=transform_matrix(k,-1,2,0,-1/637);l
array([[ 1. , 0. , -0. , 1.61538462],
[ 0. , 1. , -0. , 0.76923077],
[ 0. , 0. , 1. , -0.15384615]])
print("x1=",l[0,3]);print("x2=",l[1,3]);print("x3=",l[2,3])
x1= 1.6153846153846154
x2= 0.7692307692307689
x3= -0.15384615384615366
验证答案是否正确：
x1=1.6153846153846154
x2=0.7692307692307689
x3=-0.15384615384615366
y1=np.float64(5*x1-x2+2*x3-7)
y2=np.float64(-2*x1+6*x2+9*x3)
y3=np.float64(-7*x1+5*x2-3*x3+7)
print("y1=",y1,"y2=",y2,"y3=",y3)
y1= 0.0 y2= -6.661338147750939e-16 y3= -2.6645352591003757e-15
数组(3,4,-2)不是方程组的解 □
