Численные ошибки в расчете собственных значений (NumPy / SciPy) - PullRequest
0 голосов
/ 04 марта 2020

Я пытался получить собственные значения и собственные векторы произведения экспонент матрицы с помощью функции linalg.eig() из NumPy и столкнулся с большими числовыми ошибками.

Произошло следующее произведение:

def matrix(t): # matrix with n rows and n columns
    ...

S = scipy.linalg.expm( -1j * matrix(0) )
for i in range(N):
    S_prev = S
    S_next = scipy.linalg.expm( -1j * matrix(0 + i * delta_t) )
    S = numpy.matmul(S_next, S_prev)

Я должен заметить, что я могу не просто написать S = scipy.linalg.expm( -1j * (matrix(0) + matrix(delta_t) + ... ), потому что матрицы в показателе степени с различными t не коммутируют.

Все собственные значения S должны иметь абсолютное значение 1, поскольку в показателе степени есть мнимая единица 1j, а matrix(t) - эрмитово (все ее собственные значения действительны). Я проверяю этот факт:

vals, funcs = numpy.linalg.eig(S)
print(numpy.absolute(vals))

и вижу, что первое и последнее значения далеки от 1:

[2.01401514e+07 1.75005607e+03 1.18825323e+03 7.62941962e+02
 4.93897530e+02 3.56924391e+02 3.10300966e+02 2.54560554e+02
 1.84934167e+02 1.34093245e+02 1.25981874e+02 9.58105563e+01
 5.12598735e+01 3.87933877e+01 2.88031162e+01 1.59672649e+01
 1.18759460e+01 9.29488363e+00 8.71429012e+00 8.43350092e+00
 8.43464605e+00 4.51898108e+00 2.62851258e+00 1.71905799e+00
 1.25618921e+00 1.00000000e+00 1.00000000e+00 1.00000000e+00
 1.00000000e+00 1.00000000e+00 1.00000000e+00 1.00000000e+00
 9.88249716e-01 1.00000000e+00 1.00000000e+00 1.01188999e+00
 ...
 1.00000000e+00 1.00000000e+00 1.00000000e+00 1.00000000e+00
 9.99999997e-01 1.00000000e+00 2.21288824e-01 1.00254283e+00
 1.00036087e+00 9.99990372e-01 9.99542072e-01 1.18574662e-01
 1.18558625e-01 1.14754055e-01 1.07586057e-01 8.42037792e-02
 6.26281297e-02 3.47184828e-02 2.57776048e-02 1.95084377e-02
 7.93765396e-03 2.37035650e-08 1.04372635e-02 7.45749680e-03
 5.40732797e-03 3.92833442e-03 2.80172313e-03 3.22267950e-03
 2.02470192e-03 1.31073074e-03 8.41571767e-04 5.71394099e-04]

Также я смотрю на собственные векторы хороших собственных значений (абсолютное значение которых равно 1) и вижу совершенно неправильное поведение, такое как высокочастотные колебания:

num = numpy.arange(n)
matplotlib.pyplot.plot(num, np.absolute(funcs[200])) # good eigenvalue number is 200
matplotlib.pyplot.ylabel("Vector Component Value")
matplotlib.pyplot.xlabel("Vector Component Number")
matplotlib.pyplot.show()

График

(Я беру абсолютное значение компонентов собственного вектора, потому что оно пропорционально плотности вероятности в моей задаче)

Можно ли исправить ошибки?

Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...