Участник:Бугаков Юрий/Построение матрицы Адамара произвольного размера

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

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

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

Квадратная матрица Н порядка m с элементами ±1 называется матрицей Адамара, если выполняется условие

[math] H_m\,H_m^T = m\,E_m [/math]

Нетрудно заметить, что различные строки матрицы Адамара попарно ортогональны. Также, можно увидеть из определения, что m четно и

[math] H_{2m} = \begin{pmatrix} H_{m} & H_{m}\\ H_{m} & -H_{m}\end{pmatrix} [/math]

С матрицами Адамара связан ряд нерешенных проблем, одна из которых состоит в следующем. Мы уже видели, что порядок m матрицы Адамара при m ≥ 3 может быть лишь четным. Более того, при m ≥ 4 порядок обязан делиться на 4. До сих пор остается открытым вопрос: для любого ли m, кратного 4, существует матрица Адамара порядка m? Неизвестно, в частности, существует ли матрица Адамара порядка 268 (это наименьший порядок, кратный 4, для которого матрица Адамара еще не построена).

Часто матрицу Адамара нормализуют и рассматривают размерности степени 2 (например [2],[3] и [4]). В дальнейшем будем рассматривать только такие случаи. Переопределим матрицу Адамара следующим образом:

Матрица Адамара H_n - это матрица размера 2ⁿ × 2ⁿ, для которой справедливо равенство

[math]H_n = H_{1} \otimes H_{n-1}[/math]

[math]\begin{align} H_1 = \frac{1}{\sqrt2} &\begin{pmatrix}\begin{array}{rr} 1 & 1\\ 1 & -1 \end{array}\end{pmatrix} \end{align}[/math]

где [math] \otimes [/math] представляет собой тензорное произведение.

Представим некоторые частные примеры матриц Адамара:

[math] H_0 = +1, [/math]

[math] H_1 = \frac{1}{\sqrt2} \begin{pmatrix}\begin{array}{rr} 1 & 1\\ 1 & -1 \end{array}\end{pmatrix} [/math]

[math] H_2 = \frac{1}{2} \begin{pmatrix}\begin{array}{rrrr} 1 & 1 & 1 & 1\\ 1 & -1 & 1 & -1\\ 1 & 1 & -1 & -1\\ 1 & -1 & -1 & 1 \end{array}\end{pmatrix} [/math]

[math] H_3 = \frac{1}{2^{\frac{3}{2}}} \begin{pmatrix}\begin{array}{rrrrrrrr} 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1\\ 1 & -1 & 1 & -1 & 1 & -1 & 1 & -1\\ 1 & 1 & -1 & -1 & 1 & 1 & -1 & -1\\ 1 & -1 & -1 & 1 & 1 & -1 & -1 & 1\\ 1 & 1 & 1 & -1 & -1 & -1 & -1 & -1\\ 1 & -1 & 1 & 1 & -1 & 1 & -1 & 1\\ 1 & 1 & -1 & 1 & -1 & -1 & 1 & 1\\ 1 & -1 & -1 & -1 & -1 & 1 & 1 & -1\\ \end{array}\end{pmatrix} [/math]

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

Основой алгоритма является тот факт, что любую матрицу Адамара [math]H_n[/math], размерностью [math]2^n\times 2^n[/math] можно вычислить поэлементно ([3], [4]):

[math]h_{k,l} = \frac{1}{2^\frac{n}{2}} (-1)^{\sum_j k_j l_j}[/math],

где k_j и l_j коэффициенты двоичного разложения чисел k и l (1 или 0)

[math]k = \sum^{n}_{i=0} {k_i 2^i} = k_{n} 2^{n} + k_{n-1} 2^{n-1} + \cdots + k_1 2 + k_0[/math]

и

[math]l = \sum^{n}_{i=0} {l_i 2^i} = l_{n} 2^{n} + l_{n-1} 2^{n-1} + \cdots + l_1 2 + l_0[/math].

Помимо этой формулы чаще всего используется метод тензорного произведения. В связи с трудностями распараллеливания и постоянного выделения памяти при каждом шаге рекурсивного вызова этого метода, был выбран именно этот способ вычисления.

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

Вычислительное ядро составляют вычисления элементов матрицы [math]h_{k,l} = \frac{1}{2^\frac{n}{2}} (-1)^{\sum_j k_j l_j}[/math].

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

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

[math]\sum_{j=0}^n k_j l_j[/math].

Как видно по формуле, получаемое значение равна количеству единиц в двоичной записи результата побитового умножения k и l. Существует эффективный метод вычисления количества единиц , сложность которого равна m операций вычитания и m операций побитового умножения, где m-искомое количество единиц.

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

Последовательность исполнения метода следующая:

1. Побитовое умножение [math]k\&l[/math].

2. Вычисление количества единиц двоичной записи.

3. Определение знака элемента по степени.

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

Для построения матрицы Адамара порядка [math]N=2^n[/math] в первом шаге потребуется [math]N^2[/math] побитовых умножений k&l. Во втором шаге, как говорилось выше, сложность метода вычисления количества единиц числа в двоичном представлении равна удвоенному количеству единиц. В наихудшем случае сложность равна [math]2log_2N[/math], при k=l=N-1, а в наилучшем 0. Для нахождения степени всех элементов матрицы потребуются [math] 2N^2[\frac{log_2N}{4}][/math] операций. В третьем шаге для определения знака потребуется [math]N^2[/math] операций.

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

На рисунке ниже изображен информационный граф алгоритма.

• [math]k\&l[/math] - операция побитового умножения, где [math]k[/math] и [math]l[/math] соответствующие координаты элемента [math]h_{kl}[/math].
• OP - вычисление количества единиц числа в двоичной системе. Состоит из последовательного применения операций вычитания и побитового умножения.
• SGN - определение знака соответствующего элемента. Если результат предыдущего шага четное, то знак положительный, иначе отрицательный.

Рисунок 1. Граф алгоритма.

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

Для построения матрицы порядка [math]N=2^n[/math] требуется последовательно выполнить следующие шаги:

• [math]N^2[/math] вычислений побитового умножения,
• для каждого элемента от [math]2[/math] до [math]2n[/math] операций вычитания и побитового умножения, в зависимости от результата предыдущего шага,
• [math]N^2[/math] для определения знака значения элемента.

Таким образом для построения необходимы лишь простые операции вычитания и побитового умножения выполняемые с целыми положительными числами. Стоит заметить что вычисления для каждого элемента не зависят друг от друга. Также можно воспользоваться либо симметричностью матрицы Адамара либо закономерностью областей показанная в пункте 1.1. Таким образом, основной вклад в высотe ЯПФ вносит операция вычисления количества единиц. Ширина ЯПФ будет равна O(n^2).

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

Входные данные: n - размерность матрицы [math]H_n[/math].

Объём входных данных: 1 (число n).

Выходные данные: матрица Адамара [math]H_n[/math] размерностью [math]2^n\times 2^n[/math].

Объём выходных данных: [math]2^{2n}[/math].

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

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

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

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

Дуги информационного графа, исходящие из вершин, соответствующих операциям квадратного корня и деления, образуют пучки т. н. рассылок линейной мощности (то есть степень исхода этих вершин и мощность работы с этими данными — линейная функция от порядка матрицы и координат этих вершин). При этом естественно наличие в этих пучках «длинных» дуг. Остальные дуги локальны.

Наиболее известной является компактная укладка графа — его проекция на треугольник матрицы, который перевычисляется укладываемыми операциями. При этом «длинные» дуги можно убрать, заменив более дальнюю пересылку комбинацией нескольких ближних (к соседям).

Эквивалентное возмущение [math]M[/math] у метода Холецкого всего вдвое больше, чем возмущение [math]\delta A[/math], вносимое в матрицу при вводе чисел в компьютер: [math] ||M||_{E} \leq 2||\delta A||_{E} [/math]

Это явление обусловлено положительной определённостью матрицы. Среди всех используемых разложений матриц это наименьшее из эквивалентных возмущений.

2 ЧАСТЬ Программная реализация алгоритма

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

В простейшем варианте генерацию матрицы Адамара для нашего алгоритма можно записать так:

void adamar(int** H, int n) // N = 2^n - размерность матрицы
{     
    int i,j,N;
    float h;
    h = pow(1/2,n/2); // h - нормализующий множитель
    N = pow(2,n);

    for(i=0;i<N;++i) 
    {
	for(j=0;j<N;++j)
	{
		ij = i & j;
    		for (l=0;ij;++l) ij&=ij-1; // l - количество единиц бинарного представления ij
		if(l%2==0) H[i][j] = h;
		else H[i][j] = -h;
	}
    }
}

Также существует рекурсивная реализация метода генерации матрицы Адамара

void adamar(int** H, int N, int i, int j, int h) // h - нормализующий множитель
{
    if(n==1) H[i][j] = h;
    else
    {
        adamar(H, N/2, i, j, h);
        adamar(H, N/2, i, N/2+j, h);
        adamar(H, N/2, N/2+i, j, h);
        adamar(H, N/2, N/2+i, N/2+j, -h);
    }
}

Или эквивалентный безрекурсивный вариант:

void adamar(int** H, int n) // N = 2^n - размерность матрицы
{     
    int i,j,N;
    float h;
    h = pow(1/2,n/2);
    N = pow(2,n);
    k=2;
    H[0][0] = h;
	
    while(k<=N)
    {
	for (i=0; i<k; ++i)
	{
		for (j=0; j<k; ++j)
		{
			if (i<k/2 && j<k/2)
			{
			}
			else
			{
				if(i>=k/2 && j>=k/2)
				{
					H[i][j] = -(H[i-k/2][j-k/2]);
				}
				else
				{
					if (i>j)
					{
						H[i][j] = H[i-k/2][j];
					}
					else
					{
						H[i][j] = H[i][j-k/2];
					}
				}
			}
		}
	}
		
        k=k*2;
    }
}

2.7 Существующие реализации алгоритма

Примеры вычислительных пакетов, содержащих функции генерации матриц Адамара:

1. MATLAB
2. WOLFRAM
3. Python
4. Octave

Стоит отметить, что во многих вычислительных пакетах встроена операция тензорного произведения матриц, которая позволяет строить матрицу Адамара по определению. Этот факт позволяет использовать матрицу Адамара в этих вычислительных пакетах (MAPLE MATCATH, ...)

3 Литература

[1] Мак-Вильямс Ф.Дж., Слоэн Н.Дж.А. Теория кодов, исправляющих ошибки: Пер. с англ. — М: Связь, 1979. — 744 с.

[2] Кронберг, Ю.И. Ожигов, А.Ю. Чернявский — Алгебраический аппарат квантовой информатики 2

[3] Kunz H.O.. On the Equivalence Between One-Dimensional Discrete Walsh-Hadamard and Multidimensional Discrete Fourier Transforms. IEEE Transactions on Computers. 28 (3): 267–8. doi:10.1109/TC.1979.1675334

[4] Википедия – свободная энциклопедия [Электронный ресурс]. - https://en.wikipedia.org/wiki/Hadamard_transform. - (дата обращения: 15.10.2016).

[5] Воеводин В.В.. Воеводин Вл.В. Параллельные вычисления. - СПб.: БХВ-Петербург, 2002. - 608 с.

[6] Воронцов К.В.. Математические методы обучения по прецедентам (теория обучения машин)