主要還是調包：

from numpy.linalg import eig

特徵值分解：  A = P*B*P^T  當然也可以寫成 A = P^T*B*P  其中B為對角元為A的特徵值的對角矩陣。

首先A得正定，然後才能在實數域上分解，

>>> A = np.random.randint(-10,10,(4,4))
>>> A
array([[  6,   9, -10,  -1],
       [  5,   9,   5,  -5],
       [ -8,   7,  -4,   4],
       [ 
 
-1,  -9,   0,   6]])

>>> C = np.dot(A.T, A)
>>> C
array([[126,  52,  -3, -69],
       [ 52, 292, -73, -80],
       [ -3, -73, 141, -31],
       [-69, -80, -31,  78]])

>>> vals, vecs = eig(C)
>>> vals
array([357.33597086, 174.10172008,   8.84429957,  96.71800949])
>>> vecs
array([[
 
-0.28738314, -0.51589436, -0.38221983, -0.71075449],
       [-0.87487263,  0.12873861, -0.24968051,  0.39456798],
       [ 0.2572149 , -0.69304313, -0.33950158,  0.58161018],
       [ 0.29300052,  0.48679627, -0.82237845, -0.02955945]])

故使用時應先將特徵值轉換為矩陣：

>>> Lambda = np.diag(vals)
>>> Lambda
array([[
 
357.33597086,   0.        ,   0.        ,   0.        ],
       [  0.        , 174.10172008,   0.        ,   0.        ],
       [  0.        ,   0.        ,   8.84429957,   0.        ],
       [  0.        ,   0.        ,   0.        ,  96.71800949]])

>>> np.dot(np.dot(vecs, Lambda), vecs.T) # 與C=A.T*A相等
array([[126.,  52.,  -3., -69.],
       [ 52., 292., -73., -80.],
       [ -3., -73., 141., -31.],
       [-69., -80., -31.,  78.]])

>>> np.dot(np.dot(vecs.T, Lambda), vecs)
array([[171.65817919,  45.58778569,  53.20435074,  13.37512137],
       [ 45.58778569, 125.15670964,  28.22684299, 134.91290105],
       [ 53.20435074,  28.22684299, 129.48789571,  80.5284382 ],
       [ 13.37512137, 134.91290105,  80.5284382 , 210.69721545]])

故驗證了使用np中的eig分解為A=P*B*P^T 而不是A=P^T*B*P，其中P=vecs，

即 C = vecs * np.diag(vals) * vecs.T # 這裡簡寫*為矩陣乘法 

然後再來看使用np中的eig分解出來的vec中行向量是特徵向量還是列向量是特徵向量，只需驗證：A*vecs[0] = vals[0]*vecs[0]

>>> np.dot(C, vecs[0])
array([-12.84806258, -80.82266859,   6.66283128,  17.51094927])
>>> vals[0]*vecs[0]
array([-102.69233303, -184.34761071, -136.58089252, -253.97814676])

>>> np.dot(C, vecs[:,0])
array([-102.69233303, -312.62346098,   91.91213634,  104.69962583])
>>> vals[0]*vecs[:, 0]
array([-102.69233303, -312.62346098,   91.91213634,  104.69962583])

後者兩個是相等的，故使用np中的eig分解出的vecs的列向量是特徵向量。

然後我們可以驗證P是單位正交矩陣：

>>> np.dot(vecs.T, vecs)
array([[ 1.00000000e+00, -7.13175042e-17, -2.45525952e-18,
         2.75965773e-16],
       [-7.13175042e-17,  1.00000000e+00,  2.49530948e-17,
        -5.58839097e-16],
       [-2.45525952e-18,  2.49530948e-17,  1.00000000e+00,
        -7.85564967e-17],
       [ 2.75965773e-16, -5.58839097e-16, -7.85564967e-17,
         1.00000000e+00]])

>>> np.dot(vecs, vecs.T)
array([[ 1.00000000e+00,  2.97888865e-16, -2.68317972e-16,
         1.69020590e-16],
       [ 2.97888865e-16,  1.00000000e+00, -4.40952204e-18,
        -6.24188690e-17],
       [-2.68317972e-16, -4.40952204e-18,  1.00000000e+00,
        -1.13726775e-17],
       [ 1.69020590e-16, -6.24188690e-17, -1.13726775e-17,
         1.00000000e+00]])

# 可以看到除對角元外其他都是非常小的數

即P^T*P = P*P^T = E , P^T=P^-1。事實上，在求解P的過程中就使用了施密特正交化過程。

另一方面，我們從數學角度來看：

首先補充一些數學知識：

首先AB相似：P^-1*A*P=B，AB合同：C^T*A*C=B,

二次型：係數在K中的一個n元二次多項式。由其生成的矩陣稱為二次型的矩陣，二次型的矩陣一定是對稱矩陣！

正定矩陣：實二次型x^T*A*x > 0, x為列向量。

性質：假設A為正定矩陣

1、正定矩陣特徵值全大於0

2、行列式 |A| >0

3、A合同於單位陣E，即存在可逆方陣C， s.t. C^T*E*C = A = C^T*C, 顯然可得A為對稱正定

正交矩陣：A*A^T=A^T*A=E ，

性質：

1、A的各行/列是單位向量且兩兩正交

2、A^T=A^-1

3、|A|=1

4、(Ax,Ay)=(x,y)x,y∈R

酉矩陣：A*A^H=A^H*A=E 顯然為正交矩陣在複數域上的推廣。其中H為共軛轉置。

性質：

1、A的各行/列是單位向量且兩兩正交

2、A^H=A^-1

3、|A|=1

（這裡補充一個厄米特矩陣：A^H = A）

正規矩陣：A*A^H=A^H*A (以上的矩陣均有這個性質，故正規矩陣最為廣泛)

正規矩陣的充要條件是：存在酉矩陣U，使得A酉相似於對角矩陣B，即U^H*A*U=U^-1*A*U=B。

A = P*B*P  ，其中B為對角元素為A的特徵值的對角陣，P的列向量為特徵值對應的特徵向量(因為B每行乘以P每列)