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)


 

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

54193. ГРУПИ ТА ТИПИ ПРОСТИХ ЗАДАЧ 217 KB
  ГРУПИ ТА ТИПИ ПРОСТИХ ЗАДАЧ МЕТА: узагальнити знання студентів з теми Прості задачі вчити аналізувати добирати складати прості задачі; розвивати вміння користуватися додатковими матеріалами приміняти отриманні знання на уроках практики; розвивати логічне мислення память увагу фантазію міжпредметні звязки; виховувати самостійність швидкість реакції увагу. Тернопіль 2002р ЗАВДАННЯ НА САМОСТІЙНУ РОБОТУ Повторити конспект дібрати та проаналізувати 3 задачі з підручника Математика 3 клас на множення...
54194. АКТИВІЗАЦІЯ ПІЗНАВАЛЬНОЇ ДІЯЛЬНОСТІ УЧНІВ НА УРОКАХ МАТЕМАТИКИ З ВИКОРИСТАННЯМ ІКТ 239 KB
  Активізація пізнавальної діяльності учня без розвитку його пізнавального інтересу не тільки важка, але й практично неможлива. От чому в процесі навчання необхідно систематично збуджувати, розвивати і укріплювати пізнавальний інтерес учнів і як важливий мотив навчання, і як стійку рису особистості, і як могутній засіб виховання.
54195. Позакласний захід в 2 класі «Математичний ранок» 196.5 KB
  Мета: познайомити учнів з історією чисел, їх написанням; вчити розв’язувати завдання з логічним навантаженням, застосовуючи знання з математики, природознавства; розвивати кмітливість, спостережливість, логічне мислення учнів; виховувати інтерес до математики.
54196. НТЕЛЛЕКТУАЛЬНО-РАЗВЛЕКАТЕЛЬНАЯ ИГРА ДЛЯ УЧАЩИХСЯ 5-7 КЛАССОВ «ФАН КЛУБ МАТЕМАТИКИ» 105.5 KB
  Развитие умений формулировать и излагать мысль, моделировать ситуацию. Развитие навыков работы в группе. Воспитание стойкости, находчивости, любознательности.
54197. Клуб веселых математиков 309.5 KB
  Сегодня мы открываем клуб веселых математиков. В соревнованиях участвуют два класса - две команды. Это лучшие математики, которые не унывают, быстро считают, хорошо решают задачи, любознательны, живут весело и дружно.
54198. Сценарій позакласного заходу: «Математика – зліва, математика - справа» 134 KB
  Бажаю вам дорогі друзі шановні академіки успіхів у вивченні цариці наук Математики 1 учень вбігає Хлопці я чув що в школі буде тиждень математики Уявляєте весь тиждень сама лиш тільки математика 2 учень з місця Та не може такого бути 3 учень вбігає У мене для вас ось така новина Всіх вчителів направили на курси залишилися лише математики. 1 учень Ну. Учень Як зібрався математику вивчати То до класу на уроки вирушай. Виконується цікавий номер Учень Сорокап'ятирічний чоловік покохав п'ятнадцятирічну дівчину.
54200. Конкурс знавців математики «Мадонна Математика» 714 KB
  Сьогодні будемо ми друзі Царицю всіх наук вітати. Не всі ви в майбутньому станете математиками але математика потрібна і в науці і в техніці і в повсякденному житті. Ще в давні часи математику називали царицею наук ключем до всіх наук. Одне слово одне слово Математику Чом по курсу спішать на морях кораблі Хуртовини й тумани долають в імлі Капітани не ледачі Не лякають їх задачі Одне слово одне слово Математики Хочеш лікарем стати хочеш в космос літати Перш за все треба друже математику знати Всі повинні шанувати Ікси...
54201. Математична конференція «Золотий переріз – душа гармонії» 502 KB
  І називається вона Золотий переріз душа гармонії. Теорему Піфагора знає кожен школяр а про золотий переріз далеко не всі. Про золотий переріз знали ще в Давньому Єгипті й Вавилоні в Індії та Китаї.