我相信我的具体问题还没有得到解答。
假设我有一个 ndarray,它仅在最深的两个级别上包含需要对角化的实际方形数组。
为了简单起见,我们定义以下内容
import numpy as np
from numpy.linalg import eig
a = np.array([[1,4,8],[2,9,3],[2,6,5]])
对于这个特定矩阵,
eig(a)
给出的特征系统是
eig(a)
EigResult(eigenvalues=array([13.54300945, -1.07781094, 2.53480149]), eigenvectors=array([[-0.554087 , -0.97321006, 0.82404034],
[-0.61465955, 0.1385103 , -0.42744667],
[-0.56141004, 0.18351316, 0.37181563]]))
如您所见,特征值 13.54、-1.07 和 2.53 未排序。
现在让我们采用一个我想要对角化的实际嵌套数组系统作为测试问题
b = np.array([[a,a,a,a],[a,a,a,a]])
e,v=eig(b)
好消息,
eig()
在数组数组上进行线程化。但现在我想要两件事:首先,获取 e
上的排序,并根据该排序对 e
中的特征值和 v
中的特征向量进行排序,然后 转置特征向量,以便实际上是 v
上的最后一个轴
跨越每个特征向量的元素(现在倒数第二个轴跨越元素,我认为这是线性代数中的常见做法,其中特征向量通常以列形式给出)。
第一部分,处理
e
,可以通过组合多个解决方案来完成,即argsort
和take_along_axis
ordering=e.argsort(axis=-1)
e_sorted = np.take_along_axis(e, ordering, axis=-1)
如何对特征向量
v
进行排序,以便它们根据特征值进行排序,但同时最后两个轴被转置?无论我尝试什么,顺序都不符合 v
。
我想我已经明白了。首先,进行一些设置,以便我们有一个足够强大的示例:
import numpy as np
from numpy.linalg import eig
a = np.array([[345, 23, 2],
[23, 4, 34],
[2, 34, 5]])
b = np.array([[2, 2, 4],
[2, 3, 5],
[4, 5, 4]])
c = np.array([[a,b],[b,b]])
e,v=eig(c)
print(e[1,1])
print(v[1,1])
这里的输出是:
[10.6898736 0.43934625 -2.12921985]
[[-0.44860502 -0.77667543 -0.44218639]
[-0.56625935 0.62978926 -0.53171029]
[-0.69145056 -0.01186427 0.72232635]]
您可以验证,这是
b
的正确特征系统,但是,特征值未排序,特征向量是列(因此它们跨轴=-1)。我们想解决这两个问题。
首先,我们首先找到特征值的排序
ordering_e = e.argsort(axis=-1)
可以以非常简单的方式使用此排序来对特征值进行排序(此处
e
被修改为就地)
e[:] = np.take_along_axis(e, ordering_e, axis=-1)
对于向量,我们需要扩展排序数组的形状。具体来说,不同特征值的向量由
v
中的最后一个索引枚举,但它们的分量存储在倒数第二个索引中,我们需要在 argsort
变量中添加一个额外的虚拟维度。将其视为形状扩展:(..., dim) -> (..., dim, dim)。这可以通过多次使用 np.newaxis
来实现。次数由原始矩阵的维数给出,这也是e.shape
的最后一个维度。放在一起:
ordering_v = np.repeat(ordering_e[...,np.newaxis,:], e.shape[-1], axis=-2)
特征向量现在可以以相同的方式重新排列(也就地)
v[:] = np.take_along_axis(v, ordering_v, axis=-1)
最后可以通过numpy的
moveaxis
来完成转置。由于我们要调换最后两个级别,这样就可以了
v[:] = np.moveaxis(v, -2, -1)
注意,
moveaxis(arr, m, n)
通常与 moveaxis(arr, n, m)
不同,但由于我们要调换最后两个级别,所以 -2, -1 与 -1, -2 相同。
仅此而已,现在这两张照片给出了
[-2.12921985 0.43934625 10.6898736]
[[-0.44218639 -0.53171029 0.72232635]
[-0.77667543 0.62978926 -0.01186427]
[-0.44860502 -0.56625935 -0.69145056]]
即正是我们想要的:排序的特征值和相应排序的特征向量,其分量跨越最后一层。
TL;DR:使用此代码对特征系统进行就地
排序def sort_eigensystem(E,V):
ordering_E = E.argsort(axis=-1)
E[:] = np.take_along_axis(E, ordering_E, axis=-1)
ordering_V = np.repeat(ordering_E[...,np.newaxis,:], E.shape[-1], axis=-2)
V[:] = np.moveaxis(np.take_along_axis(V, ordering_v, axis=-1), -2, -1)