Однокубитное преобразование вектора-состояния

Однокубитное преобразование вектора-состояния
Последовательный алгоритм
Последовательная сложность	[math]3 \cdot 2^n[/math]
Объём входных данных	[math]2^n+4[/math]
Объём выходных данных	[math]2^n[/math]
Параллельный алгоритм
Высота ярусно-параллельной формы	[math]2[/math]
Ширина ярусно-параллельной формы	[math]2^{n+1}[/math]

Основные авторы описания: А.Ю.Чернявский.

1 Свойства и структура алгоритма

1.1 Общее описание алгоритма

Алгоритм производит моделирование действия однокубитного квантового вентиля на вектор-состояние. Данный алгоритм обычно является подпрограммой и многократно применяется к различным кубитам одного состояния (например при моделировании квантовых алгоритмов или анализе квантовой запутанности). Особенностью алгоритма, как и большинства алгоритмов квантовой инфоорматики, является экспоненциальный рост объема данных в зависимости от основного параметра - числа кубитов, что приводит к необходимости суперкомпьютерной реализации для решения важных практических задач.

1.2 Математическое описание алгоритма

Исходные данные:

Целочисленные параметры [math]n - [/math] число кубитов (необязательно) и [math]k -[/math] номер кубита, над которым производится преобразование.

Комплексная матрица [math]U = \begin{pmatrix} u_{00} & u_{01}\\ u_{10} & u_{11} \end{pmatrix}[/math] однокубитного преобразования размера [math]2 \times 2.[/math]

Комплексный вектор [math]v[/math] размерности [math]2^n,[/math] задающий начальное состояние многокубитной системы.

Вычисляемые данные: комплексный вектор [math]w[/math] размерности [math]2^n,[/math] соответствующий состоянию после преобразования.

Формулы метода:

Состояние после действия преобразования [math]U[/math] на [math]k-[/math]й кубит имеет вид [math]v_{out} = I_{2^{k-1}}\otimes U \otimes I_{2^{n-k}},[/math] где [math]I_{j} - [/math] единичная матрица размерности [math]j,[/math] а [math]\otimes - [/math] тензорное произведение (произведение Кронекера).

Однако, элементы итогового вектора можно записать и в прямом виде, что более удобно для вычислений:

[math] w_{i_1i_2\ldots i_k \ldots i_n} = \sum\limits_{j_k=0}^1 u_{i_k j} v_{i_1i_2\ldots j_k \ldots i_n} = u_{i_k 0} v_{i_1i_2\ldots 0_k \ldots i_n} + u_{i_k 1} v_{i_1i_2\ldots 1_k \ldots i_n} [/math]

Индекс-кортеж [math]i_1i_2\ldots i_n[/math] представляет собой двоичную запись индекса элемента в массиве.

1.3 Вычислительное ядро алгоритма

Вычислительное ядро алгоритма представляет собой независимое вычисление всех [math]2^n[/math] элементов вектора [math]w.[/math] Вычисление каждого элемента требует две операции умножения и одну операцию сложения. Кроме того необходимо вычислять индексы типа [math]i_1i_2\ldots 0_k \ldots i_n,[/math] а также значение бита [math]i_k,[/math] что требует побитовых операций.

1.4 Макроструктура алгоритма

Как записано и в описании ядра алгоритма, основную часть метода составляют независимые вычсиления элементов выходного вектора.

1.5 Схема реализации последовательного алгоритма

Для индекса [math]i[/math] от [math]0[/math] до [math]2^n-1[/math]

1. Вычислить элемент [math]i_k[/math] двоичного представления индекса [math]i.[/math]
2. Вычислить индексы [math]j[/math] имеющие двоичные представления [math]i_1i_2\ldots \overline{i_k} \ldots i_n,[/math] где крышка означает обращение бита.
3. Вычислить [math]w_i = u_{i_k i_k}\cdot v_{i} + u_{i_k \overline{i_k}}\cdot v_j.[/math]

1.6 Последовательная сложность алгоритма

Алгоритм требует:

1. [math]2^{n+1}[/math] операций умножения комплексных чисел;
2. [math]2^n[/math] операций сложения комплексных чисел;
3. [math]2^n[/math] операций получения значения [math]k[/math]-го бита числа;
4. [math]2^n[/math] операций изменения значения [math]k[/math]-го бита числа.

Отметим, что данный алгоритм обычно применяется много раз подряд, в связи с чем вычисления, связанные с побитовыми операциями (3-4), могут однократно проводиться в начале алгоритма. Кроме того, от них можно избавиться, пользуясь сложением и логическим умножением с числом [math]2^k,[/math] которое сохраняется для всего алгоритма.

1.7 Информационный граф

Представим рисунки графов алгоритма для случая [math]n=3, k=1[/math] (рис.1) и [math]n=3, k=2[/math] (рис.2). На графах не представлены матрицы преобразования [math]U,[/math] в связи с тем, что их размер при больших [math]n[/math] много меньше, нежели размеры входного и выходного векторов. На Рис.3 изображена основная операция, представляемая на Рис.1 и Рис.2 оранжевым цветом.

Отображение графа со входными и выходными данных, а также "свёрнутой" тройной операцией удобно для понимания локальности обращений к памяти. Отметим, что структура графа (а именно обращение к входным данным) сильно зависит от параметра [math]k.[/math]

Рисунок 1. Граф алгоритма для [math]n=3, k=1[/math] без отображения матрицы преобразования [math]U.[/math]

Рисунок 2. Граф алгоритма для [math]n=3, k=2[/math] без отображения матрицы преобразования [math]U.[/math]

Рисунок 3. Граф основной операции.

1.8 Ресурс параллелизма алгоритма

Как видно из информационного графа, прямой алгоритм моделирования однокубитного преобразования обладает высочайшей степенью параллелизма. Все операции вычисления элементов нового вектора-состояния могут быть произведены параллельно. Для вычисления одного элемента необходимо выполнить две операции умножения и одну операцию сложения, операции умножения, в свою очередь, также могут быть выполнены параллельно. Таким образом, необходимо выполнение двух ярусов: одного, состоящего из [math]2^{n+1}[/math] умножений и другого, состоящего из [math]2^n[/math] сложений.

1.9 Входные и выходные данные алгоритма

Входные данные:

• Комплексный вектор состояния [math]u[/math] длины [math]2^n.[/math] Обычно нормирован на единицу.
• Унитарная матрица [math]U[/math] порядка [math]2[/math].
• Номер кубита [math]q[/math].

Выходные данные:

• Комплексный вектор состояния [math]w[/math] длины [math]2^n[/math].

Объем входных данных: [math]2^n+4[/math] комплексных чисел и [math]1[/math] целочисленный параметр.

Объем выходных данных: [math]2^n[/math] комплексных чисел.

1.10 Свойства алгоритма

Соотношение последовательной и параллельной сложности является экспоненциальным (эксмпонента переходит в константу).

Вычислительная мощность алгоритма, как отношение числа операций к суммарному объему входных и выходных данных – константа.

2 Программная реализация алгоритма

2.1 Особенности реализации последовательного алгоритма

На языке C функцию однокубитного преобразования можно записать следующим образом:

void OneQubitEvolution(complexd *in, complexd *out, complexd U[2][2], int nqubits, int q)
{
	//n - число кубитов
        //q - номер кубита для преобразования

	int shift = nqubits-q;
	//Все биты нулевые, кроме соответствующего позиции преобразуемого кубита  
	int pow2q=1<<(shift);

	int N=1<<nqubits;
	for	(int i=0; i<N; i++)
	{
		//Обнуления меняющегося бита
		int i0 = i & ~pow2q;
		
		//Установка меняющегося бита
		int i1 = i | pow2q;
		
		//Получение значения меняющегося бита
		int iq = (i & pow2q) >> shift;

		out[i] = U[iq][0] * in[i0] + U[iq][1] *in[i1];
	}
}

Отметим, что существенная часть вычислений и логики кода приходится на битовые операции, однако, этого можно избежать: однокубитное преобразование в большинстве случаев является лишь подпрограммой и применяется к разным кубитам большое число раз. В свою очередь, вычисляемые при помощи битовых операций индексы i0, i1 и iq зависят лишь от параметров количества кубитов n и номера кубита q. Число кубитов обычно фиксировано, соответственно, можно вычислить эти индексы заранее для всех q от 1 до n. Для хранения потребуется лишь массив целочисленных переменных линейного размера 3n в то время, как обрабатываемые данные имеют экспоненциальный размер. Очевидно, что такая оптимизация критически необходима при реализации алгоритма на вычислительных устройствах или языках программирования с отсутствием быстрых битовых операций (примером может служить среда Matlab).

2.2 Возможные способы и особенности параллельной реализации алгоритма

Основной способ параллельной реализации очевиден - необходимо распараллеливание основного цикла (параллельное вычисление различных компонент выходного вектора-состояния) и, желательно, операций умножения. На машинах с общей памятью такой вариант распараллеливания приводит к ускорению, близкому к максимально-возможному. Однако, данный способ сталкивается с проблемами локальности данных. Анализируя математическая описание и информационные графы алгоритма легко видеть, что при использовании большого числа узлов с собственной памятью количество необходимых пересылок между узлами становится сопоставимым с количеством вычислений, что приводит к существенной потере эффективности. Возможны разные пути частичного решения этой проблемы, например, кэширование или использование парадигмы программирования SHMEM, однако, столь сильная нелокальность использования данных всё равно не позволяет добиться хорошей эффективности.

Алгоритм реализован (в основном последовательные версии) в различных библиотеках для квантовой информатики и квантового компьютерного моделирования. Например: QLib, libquantum, QuantumPlayground, LIQUiD.

2.4 Выводы для классов архитектур

Исходя из высочайшей возможности параллелизации и при этом наличия существенной нелокальности обращения к данным, эффективная и хорошо масштабируемая параллельная реализация алгоритма однокубитного квантового преобразования легко достижима на машинах с общей памятью. Реализация же на машинах с разделяемой памятью имеет низкую эффективность и требует специальных подходов для уменьшения времени работы. Можно отметить, что данная задача является хорошим плацдармом для разработки методов решения задач с интенсивным использованием данных и низкой локальностью.