Быстрые собственные векторы нахождения на КВД с использованием Qutip и slepc4py

Question

Я делаю некоторые числовые моделирования квантовых вычислений, и я хотел бы найти собственные векторы большой эрмитовой матрицы (~ 2^14 строк/столбцов)

Я нахожусь работающих на 24-ядерном/48-канальном потоке XEON. Код был первоначально написан с помощью библиотеки Qutip. Я узнал, что включенная функция eigenstates() использует только один поток на моей машине, поэтому я пытаюсь найти более быстрый способ сделать это.

Я попытался с помощью scipy.linalgeig() и eigh() функции, а также scipy.sparse.linalgeig() и eigh(), но оба, кажется медленнее, чем функция, построенной в Qutip.

Я видел некоторое предположение, что я могу получить некоторое ускорение от использования файла slepc4py, однако документация пакета кажется очень недостаточной. Я не могу понять, как преобразовать комплексный массив numpy в матрицу SLEPC.

A = PETSc.Mat().create() 
A[:,:] = B[:,:] 
# where B is a scipy array of complex type 

TypeError: Cannot cast array data from dtype('complex128') to dtype('float64') according to the rule 'safe' 

Answer 1

eigensolver в QuTiP использует SciPy eigensolver. Сколько потоков используется, зависит от библиотеки BLAS, к которой привязан SciPy, а также от того, используете ли вы разреженный или плотный решатель. В плотном случае eigensolver будет использовать несколько ядер, если базовый BLAS будет использовать преимущества (например, Intel MKL). Редкий разрешитель использует в основном разреженные операции matvec, которые ограничены пропускной способностью памяти и, таким образом, наиболее эффективны с использованием одного ядра. Если вы хотите все собственные значения, то вы в основном застреваете с использованием плотных решателей. Однако, если вам нужно только несколько. Такие, как самые низкие несколько собственных состояний, а затем разреженные - это путь.

Answer 2

В итоге я нашел более простой способ использования всех ядер, похоже, что qutip не сказал mkl использовать все ядра. в моем питоне коде, я добавил:

import ctypes 
mkl_rt = ctypes.CDLL('libmkl_rt.so') 
mkl_get_max_threads = mkl_rt.mkl_get_max_threads 
mkl_rt.mkl_set_num_threads(ctypes.byref(ctypes.c_int(48))) 

это заставил Intel MKL использовать все ядра, и дал мне хорошее ускорение.

(ответ от question)