39413

Реализация и исследование быстрого алгоритма двумерного вещественного ДПФ по основанию 4 представлением данных в гиперкомплексной алгебре

Курсовая

Информатика, кибернетика и программирование

Заданный алгоритм был реализован программно с помощью технологии Microsoft. NET Framework на языке программирования C++. Написанное приложение состоит из двух сборок: библиотеки классов FFT, содержащей все необходимое для вычисления ДПФ по формуле и БПФ.

Русский

2013-10-04

294.73 KB

20 чел.

Федеральное агентство по образованию

Государственное образовательное учреждение

высшего профессионального образования

«Самарский государственный аэрокосмический

университет  имени академика С.П. Королева»

Кафедра геоинформатики

ПОЯСНИТЕЛЬНАЯ ЗАПИСКА

к курсовому проекту по дисциплине

"Компьютерная алгебра"

Вариант № 11

Cтудент: Терентьева Д.А.
группа 65
7

Руководитель проекта: Чичева М. А.

Самара 2010

Задание

Реализация  и исследование быстрого алгоритма двумерного вещественного ДПФ по основанию 4 представлением данных в гиперкомплексной алгебре.


Содержание

1 Постановка задачи 4

2 Теоретический материал 5

2.1 Двумерное преобразование Фурье 5

2.1.1 Определение двумерного ДПФ в алгебре комплексных чисел  и в четырехмерных алгебрах 5

2.1.2 Алгоритм двумерного ДПФ по основанию 2 6

2.1.3 Алгоритм двумерного ДПФ по основанию 4 7

2.2 Алгебра гиперкомплексных чисел 7

2.3 Прямая сумма комплексных алгебр 9

3 Описание разработанного алгоритма 10

3.1 Описание работы программы 10

4 Исследование алгоритма 12

Приложение А. Текст программы 13

Приложение Б. Принцип работы программы 19


1 Постановка задачи

Требуется разработать и реализовать быстрый алгоритм двумерного преобразования Фурье в гиперкомплексной алгебре по основанию 4.

На входе – квадратная матрица размера , где где .

Тестирование полученной реализации алгоритма, ее исследование и сравнение с алгоритмом двумерного ДПФ.


2 Теоретический материал

2.1 Двумерное преобразование Фурье

2.1.1 Определение двумерного ДПФ в алгебре комплексных чисел
и в четырехмерных алгебрах

Пусть - входной массив размера N×N, тогда его комплексный спектр Фурье

 ,  , (2.1)

где - комплексный корень из единицы степени N.

Под дискретным преобразованием Фурье со специальным представлением данных в рамках курсового проекта будем понимать преобразование вида

 . (2.2)

Основная идея такого преобразования заключается в погружении входного сигнала и корней в четырехмерную алгебру (в рамках курсового проекта это одна из алгебраических структур, описанных в разделе 1). При этом корни лежат в разных экземплярах поля комплексных чисел, вложенных в , у них разные мнимые единицы - i и j.

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

Поскольку элемент алгебры (кватернион или гиперкомплексное число ) определяется набором четырех вещественных чисел , то комплексный спектр (2.1) может быть получен из спектра (2.2) в алгебре   следующим образом:

 

где - вектор компонент спектра,

 .

Таким образом, мультипликативная сложность вычисления совпадает со сложностью вычисления спектра в алгебре , т.к. умножения на матрицы L, I не требуют выполнения нетривиальных операций вещественного умножения.

Далее излагаются схемы декомпозиции двумерного дискретного преобразования Фурье (2.2) и способы его вычисления.

2.1.2 Алгоритм двумерного ДПФ по основанию 2

Пусть . Представим (2.2) в виде четырех сумм, разделяя входную последовательность по четным и нечетным значениям каждого индекса :

  (2.3)

где

Вычисление спектра для остальных значений пар производится без дополнительных умножений и может быть записано в матричной форме

 

или в следующем виде

 , (2.4)
.


2.1.3 Алгоритм двумерного ДПФ по основанию 4

Разобьем теперь входную последовательность на 16 подпоследовательностей и запишем (2.2) в виде

  (2.5)

где

 , (2.6)

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

 , (2.7)
.

Умножения на степени базовых элементов i и j тривиальны, они сводятся к перестановкам элементов кода и/или смене знака компонент.

2.2 Алгебра гиперкомплексных чисел

Коммутативной гиперкомплексной алгеброй будем называть четырехмерную ассоциативную алгебру над

  (2.2.11)

с определяющими соотношениями для умножений базисных элементов :

  (2.2.12)

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

  (2.2.13)

Отметим еще два полезных свойства алгебры :

- поле комплексных чисел  канонически вкладывается в :

 


- справедливы соотношения

 , (2.2.14)
,  .

Непосредственное умножение элементов алгебры по формуле (2.2.13) требует 16 вещественных умножений и 12 вещественных сложений. Как обычно, будем считать, что:

- при оценке вычислительной сложности один из сомножителей предполагается постоянным и, следовательно, все арифметические операции над его компонентами могут быть реализованы заранее;

- умножения на степени числа 2 являются более простыми операциями, чем сложения и умножения, и не учитываются при анализе вычислительной сложности алгоритмов ДПФ.

Пусть - элемент алгебры общего вида; - комплексное число, являющееся константой в контексте рассматриваемых алгоритмов. В соответствии с представлением (2.2.14) вычисление произведения эквивалентно выполнению двух комплексных умножений. При использовании схемы "три умножения, три сложения" искомое произведение принимает вид

  (2.2.15)

и вычисляется посредством шести вещественных умножений и шести вещественных сложений.

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

Непосредственной проверкой легко убедиться также в справедливости следующего утверждения.

Отображения

  (2.2.16)

сохраняют сумму и произведение элементов алгебры , действуют тождественно на (то есть являются автоморфизмами над ).


2.3 Прямая сумма комплексных алгебр

Известно, что коммутативно-ассоциативная гиперкомплексная алгебра (2.2.11)-(2.2.12) изоморфна прямой сумме комплексных алгебр:

 .

Для перехода от стандартного представления гиперкомплексного числа

  (2.3.1)

к его представлению в прямой сумме комплексных алгебр, выполним замену переменных:

 , , , . (2.3.2)

В этом случае произвольный элемент алгебры примет вид

 

Таблица умножения новых базисных элементов имеет вид, показанный ниже.

Таблица 1 Правила умножения новых базисных элементов

Базисные
элементы

0

0

0

0

0

0

0

0

В таком представлении произведение двух произвольных элементов алгебры

 

разбивается на два независимых произведения, подобных произведениям комплексных чисел

  ,

и требует 6 вещественных умножений и 6 вещественных сложений.

Все свойства гиперкомплексной алгебры в таком представлении сохраняются.


3 Описание разработанного алгоритма

3.1 Описание работы программы

Заданный алгоритм был реализован программно с помощью технологии Microsoft .NET Framework на языке программирования C++.

Написанное приложение состоит из двух сборок: библиотеки классов FFT, содержащей все необходимое для вычисления ДПФ по формуле и БПФ.

Также присутствует консольный модуль, реализующий пользовательский интерфейс:

Рисунок 1 Основное меню программы

После выбора алгоритма пользователь видит консоль с предложением ввести размер квадратного массива:

Рисунок 2 Консоль ввода массива

В результате генерируется файл с подаваемым на вход сигналом, на выходе получаем файл с исходным спектром в комплексном представлении. Далее представлены примеры выполненной работы алгоритмов БПФ и ДПФ для длины преобразования 64:

Рисунок 3 Пример работы алгоритма вычисления БПФ-64

Рисунок 4 Пример работы алгоритма вычисления ДПФ-64


4 Исследование алгоритма

Для наглядного представления приведем таблицу и график сравнения алгоритмов БПФ и ДПФ.

Проведем сравнение результатов непосредственного вычисления дискретного преобразования Фурье и результатов быстрых алгоритмов для длин преобразований 2, 4, 8, 16, 32, 64, 128, 256.

Таблица 2 Сравнение времени выполнения БПФ и  ДПФ

N

Время вычисления БПФ, с

Время вычисления ДПФ, с

8

0,001322

0,001686

16

0,007954

0,018803

32

0,047996

0,455932

64

0,287401

8,510700

128

1,714250

160,6980

256

8,678940

3559,830

 

Рисунок 5 Зависимость времени выполнения алгоритма от длины преобразования

Из Рисунка 5 видно, что быстрый алгоритм Фурье имеет значительное преимущество в скорости перед непосредственным дискретным преобразованием, хотя для небольшого объёма данных различия несущественны. Так же присутствует погрешность, связанная с машинными характеристиками, которая напрямую влияет на скорость алгоритмов.


Приложение А. Текст программы

Class Hypercomplex

#include "StdAfx.h"

#include "Hypercomplex.h"

#include <iostream>

#include <math.h>

#include <float.h>

Hypercomplex::Hypercomplex()

{

 W = 0.0;//1.0

 X = 0.0;

 Y = 0.0;

 Z = 0.0;

}

Hypercomplex::Hypercomplex(const double w, const double x, const double y, const double z)

{

 W = w;

 X = x;

 Y = y;

 Z = z;

}

void Hypercomplex::PrintOn(ostream &strm) const

{

Hypercomplex q(W, X, Y, Z);

 if(abs(q.getW()) < 0.000001)

{

 q.setW(0.0);

}

 if(abs(q.getX()) < 0.000001)

{

 q.setX(0.0);

}

 if(abs(q.getY()) < 0.000001)

{

 q.setY(0.0);

}

 if(abs(q.getZ()) < 0.000001)

{

 q.setZ(0.0);

}

strm << "(" << q.W << "," << q.X << "," << q.Y << "," << q.Z << ")";

}

Hypercomplex operator * (const Hypercomplex &a, const Hypercomplex &b)

{

 double w,x,y,z;

 w = a.W*b.W - (a.X*b.X + a.Y*b.Y + a.Z*b.Z);

 

 x = a.W*b.X + b.W*a.X + a.Y*b.Z + a.Z*b.Y;

 y = a.W*b.Y + b.W*a.Y + a.Z*b.X + a.X*b.Z;

 z = a.W*b.Z + b.W*a.Z + a.X*b.Y + a.Y*b.X;

 return Hypercomplex(w,x,y,z);

}

const Hypercomplex& Hypercomplex::operator *= (const Hypercomplex &q)

{

 double w = W, x = X, y = Y, z = Z;

 double qW = q.W, qX = q.X, qY = q.Y, qZ = q.Z;

W = w*qW - (x*qX + y*qY + z*qZ);

X = w*qX + qW*x + y*qZ + z*qY;

Y = w*qY + qW*y + z*qX + x*qZ;

Z = w*qZ + qW*z + x*qY + y*qX;

 return *this;

}

Hypercomplex operator + (const Hypercomplex &a, const Hypercomplex &b)

{

 double w,x,y,z;

w = a.W+b.W;

x = a.X+b.X;

y = a.Y+b.Y;

z = a.Z+b.Z;

 return Hypercomplex(w,x,y,z);

}

const Hypercomplex& Hypercomplex::operator += (const Hypercomplex &q)

{

W = W+q.W;

X = X+q.X;

Y = Y+q.Y;

Z = Z+q.Z;

 return *this;

}

Hypercomplex operator - (const Hypercomplex &a, const Hypercomplex &b)

{

 double w,x,y,z;

w = a.W-b.W;

x = a.X-b.X;

y = a.Y-b.Y;

z = a.Z-b.Z;

 return Hypercomplex(w,x,y,z);

}

const Hypercomplex& Hypercomplex::operator -= (const Hypercomplex &q)

{

W = W-q.W;

X = X-q.X;

Y = Y-q.Y;

Z = Z-q.Z;

 return *this;

}

Hypercomplex Hypercomplex::operator = (const Hypercomplex &q)

{

W = q.W;

X = q.X;

Y = q.Y;

Z = q.Z;

 return *this;

}

const Hypercomplex& Hypercomplex::operator ~ ()

{

X = -X;

Y = -Y;

Z = -Z;

 return *this;

}

const Hypercomplex& Hypercomplex::exp()

{

 double Mul;

 double Length = sqrt(X*X + Y*Y + Z*Z);

 if (Length > 1.0e-4)

Mul = sin(Length)/Length;

 else

Mul = 1.0;

W = cos(Length);

X *= Mul;

Y *= Mul;

Z *= Mul;

 return *this;

}

const double& Hypercomplex::arg()

{

 return acos(W/sqrt(W*W + X*X + Y*Y + Z*Z));

}

const Hypercomplex& Hypercomplex::sgn()

{

Hypercomplex q = *this;

 return (1/sqrt(W*W + X*X + Y*Y + Z*Z))*q;

}

const Hypercomplex& Hypercomplex::log()

{

 double Length;

 Length = sqrt(X*X + Y*Y + Z*Z);

 Length = atan(Length/W);

 W = 0.0;

 X *= Length;

 Y *= Length;

 Z *= Length;

 return *this;

}

const Hypercomplex& Hypercomplex::u(Hypercomplex& q)

{

Hypercomplex u = q;

u.W = 0;

 return u;

}

const Hypercomplex Hypercomplex::pow(int k)

{

Hypercomplex b(1,0,0,0);

Hypercomplex a = *this;

 while (k)

{

 if (k & 1)

 {

  b *= a;

 }

 a *= a;

 k >>= 1;

}

 return b;

}

const Hypercomplex Hypercomplex::toHypercomplex(const double &d)

{

 double w, x, y, z;

w = d;

x = 0;

y = 0;

z = 0;

 return Hypercomplex(w,x,y,z);

}

const Hypercomplex Hypercomplex::toHypercomplex(const std::complex<double> &c)

{

 double w,x,y,z;

w = c.real();

x = c.imag();

y = 0;

z = 0;

 return Hypercomplex(w,x,y,z);

}

const Hypercomplex& Hypercomplex::Initialize(const double& w, const double& x, const double& y, const double& z)

{

W = w;

X = x;

Y = y;

Z = z;

 return *this;

}

void Hypercomplex::setW(const double &value)

{

W = value;

}

void Hypercomplex::setX(const double &value)

{

X = value;

}

void Hypercomplex::setY(const double &value)

{

Y = value;

}

void Hypercomplex::setZ(const double &value)

{

Z = value;

}

double Hypercomplex::getW()

{

 return W;

}

double Hypercomplex::getX()

{

 return X;

}

double Hypercomplex::getY()

{

 return Y;

}

double Hypercomplex::getZ()

{

 return Z;

}

Class DPF

#pragma once

#include <iostream>

#include <fstream>

#include <complex>

#include "Hypercomplex.h"

#include <conio.h>

#include <math.h>

using namespace std;

#ifndef PI

#define PI 3.14159265358979323846

#endif

class DPF

{

private:

 int N1, N2;

Hypercomplex w1, w2;

Hypercomplex** mas_w;

public:

DPF();

~DPF(void);

 double** loadMatrix(char name[]);

 void writeComplexMatrix(char name[], complex<double>** complex_matrix, int&);

 

 template<typename T>

T* binInverse(T* x, int& n)

{

 int dl = n/2;

 int st = n - 1;

 int j = 0;

 T s0 = 0;

 int k = 0;

 for(int i = 0; i < st; i++)

 {

  if(i < j)

  {

   s0 = x[i];

   x[i] = x[j];

   x[j] = s0;

  }

  k = dl;

  while(k <= j)

  {

   j -= k;

   k >>= 1;

  }

  j += k;

 }

 return x;

}

 template<typename T>

T** binMtrInverse(T** matr, int& N)

{

 int size = 0;

  size = N;

  for(int i = 0; i < size; i++)

  {

   binInverse(matr[i], size);

  }

  T* v = new T[size];

  for(int i=0; i<size; i++)

  {

   for(int j=0; j<size; j++)

   {

    v[j] = matr[j][i];

   }

   binInverse(v, size);

   for(int j=0; j<size; j++)

   {

    matr[j][i] = v[j];

   }

  }

  delete v;

 return matr;

}

Hypercomplex** matrInv(Hypercomplex** q, int& N);

 void genMatr(int N, char name[]);

Hypercomplex** _dpf(Hypercomplex** matrix);

Hypercomplex** Q44(Hypercomplex&, Hypercomplex&,  Hypercomplex&, Hypercomplex&, Hypercomplex&, Hypercomplex&,  

Hypercomplex&, Hypercomplex&,  Hypercomplex&, Hypercomplex&,  Hypercomplex&, Hypercomplex&,  Hypercomplex&, Hypercomplex&,  Hypercomplex&, Hypercomplex&);

Hypercomplex butterfly2(Hypercomplex&, Hypercomplex&, Hypercomplex&, Hypercomplex&, int&, int&, int&);

Hypercomplex Q_dft_basis4(Hypercomplex**, int&, const int&, const int&);

 int pos(int a);

Hypercomplex** HyperMatrix(const int&);

Hypercomplex** toHyperMatrix(double**, int&);

Hypercomplex** Q_ftfbs(Hypercomplex**, int &size);

Hypercomplex Q00(Hypercomplex**, int&, int&, int&);

Hypercomplex Qab(Hypercomplex**, int&, int&, int&, const int&, const int&);

 void shift(Hypercomplex**, int, int);

 void clcCoordinates(int&, int&, int, int);

 void clcMasW(int&);

 void clearMasW();

Hypercomplex** Q_algftf(Hypercomplex** , int&);

 

 void fft2(Hypercomplex**,int&);

 void fft4(Hypercomplex**,int&);

complex<double> toComplex(Hypercomplex&);

 const Hypercomplex automorphism(Hypercomplex&, const int&);

complex<double>** allocPartSpectr(Hypercomplex**, int&);

complex<double>** getComplexSpectr(Hypercomplex**, int&);

 void set(int& N);

 int get_N1();

 int get_N2();

 void setW1(int);

 void setW2(int);

};


Приложение Б. Принцип работы программы

На вход подаётся файл с числами, сгенерированный таким образом:

void DPF::genMatr(int N, char name[])

{

ofstream out(name);

out << N << ' ' << N << endl;

for(int k = 0; k < N; k ++ )

{

 for(int j = 0; j < N; j ++ )

 {

  out << pow(3,double(k*j % 3)) << ' ';

 }

 out << endl;

}

 

out.close();

}.

Рассмотрим входной сигнал размером 4х4. Текстовый файл signal.txt имеет структуру:

4 4

1 1 1 1

1 3 9 1

1 9 3 1

1 1 1 1

Получаем результат в комплексной алгебре.

Результат БПФ по основанию 4, записанный в файл spectr_bpf.txt

(6,0) (10,2) (-2,0) (10,-2)

(10,2) (2,2) (-8,0) (4,0)

(-2,0) (-8,0) (2,0) (-8,0)

(10,-2) (4,0) (-8,0) (2,-2)

Результат ДПФ по основанию 4, записанный в файл spectr_dpf.txt

(36,0) (-10,10) (0,0) (-10,-10)

(-10,10) (0,-16) (6,6) (4,0)

(0,0) (6,6) (-12,0) (6,-6)

(-10,-10) (4,0) (6,-6) (0,16)


 

А также другие работы, которые могут Вас заинтересовать

84337. Розробка зв’язного радіоприймача УКХ діапазону 808.89 KB
  В результаті проектування були розкриті наступні питання, які стосуються перспективи розвитку курсового проекту. При вивченні методів побудови аналогів проектованого пристрою і порівняння їх переваг та недоліків, були сформовані технічні умови для радіоприймача КХ діапазону
84338. Основы окружающей среды и энергосбережения 37.39 KB
  Раскройте сущность административного управления природопользованием и охраной окружающей среды. Изложите основные направления Национальной стратегии устойчивого развития Республики Беларусь в области охраны окружающей среды и энергосбережения. В процессе решения вопроса эффективного и безопасного для окружающей среды применения пестицидов реализуются разные подходы.
84339. Ассортимент и качество рыбных консервов в масле 51.21 KB
  Рыбные консервы в масле являются полноценным пищевым продуктом имеющим оригинальную рецептуру приготовления. Рыбные консервы в заливке на основе масла составляют классический ассортимент где проявляются все оттенки вкуса и запаха рыбы. Именно благодаря этим свойствам консервы являются одним...
84340. Работа с электронными геодезическими приборами 1.18 MB
  Прибор удобен в обращении благодаря быстроте производства угловых и линейных измерений широкому диапазону измеряемых расстояний возможности производить измерения без использования отражателя удобному настраиваемому под требования пользователя меню.
84341. Характеристика продукции сетевой экономики. Интернет как инструмент маркетинга 167.5 KB
  Изобретение телеграфа телефона радио и компьютера привело к беспрецедентной интеграции всех этих средств в глобальную информационную гипермедийную систему которая одновременно является средой для сотрудничества и общения людей средством общемирового вещания и распространения информации а также мощным инструментом для ведения бизнеса лишенного...
84342. Техника и технология вызова притока заменой бурового раствора (промывки) 644 KB
  После выхода из турбобура буровой раствор омывает забой скважины и как и при роторном способе бурения выносит продукты бурения на устье скважины. Сложность этой конструкции состоит в том что буровой раствор на забой скважины должен проходить через электродвигатель и редуктор. Недостатками турбобуров являются высокая чувствительность к вязкости бурового раствора и высокая частота вращения которая приводит к повышенной разработке ствола скважины при бурении мягких пород а также ускоренному износу ПРИ и следовательно к увеличению...
84343. Поводы и основания к возбуждению уголовного дела. Процессуальный порядок возбуждения и отказа в возбуждении уголовного дела 134.58 KB
  Поводы и основания к возбуждению уголовного дела. Процессуальный порядок возбуждения и отказа в возбуждении уголовного дела Поводами для возбуждения уголовного дела служат: ст. Основанием для возбуждения уголовного дела является наличие достаточных данных указывающих на признаки преступления. Заявление о преступлении может быть сделано в устном или письменном виде.
84344. Финансовое регулирование в системе управления финансами 65.89 KB
  Финансовое регулирование осуществляется в целях корректировки со стороны государства развития общественного производства в нужном направлении. Таким образом финансовое регулирование представляет собой организованную государством деятельность по использованию всех аспектов финансовых отношений в целях...
84345. Функции культуры 98 KB
  Характерно что проблематика различных функций культуры разрабатывалась разными учеными сосредоточившими свое концептуальное внимание на тех или иных сторонах духовной системы. Можно выделить следующие функции культуры: а поддержание преемственности традиции...