Алговики - Вклад [ru]

Участник:Kiseliov/Метод регуляризации Тихонова

2022-11-03T10:38:06Z

Kiseliov: /* Существующие реализации алгоритма */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительным ядром являются вычисления по формуле из пункта 1.2 в позиции <math>(i, j)</math>, поскольку формула применяется к позиции <math>(i, j)</math> независимо от других точек.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на псевдокоде в стиле C:

<syntaxhighlight lang="C">
tikhonov_regularization(N, *image, *kernel, R, alpha, r, *restored)
{
// ДельтаX в формуле (логично посчитать эти значения единожды)
delta_x = (4 * R) / N;
delta_x_pow2 = delta_x * delta_x;
delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Это выражение - константа, и его логично посчитать заранее. CONSTANT_PI - число пи
constant = (CONSTANT_PI / (2 * R)) * (CONSTANT_PI / (2 * R));

// Цикл по позициям (i, j)
for (i = 0; i < N; i++)
for (j = 0; j < N; j++)
{
// Предполагается, что матрицы хранятся в виде 1D массива, поэтому нужно вычислить правильный индекс
index = i * N + j;

// Мы сдвигаем центр системы координат в центр изображения
m = i - N / 2;
n = N / 2 - j;

// Вычисление формулы. Предполагается, что pow возводит число в степень, conj возвращает сопряжённое комплексное число,
// .real() и .imag() - получение действительной или мнимой частей комплексного числа соответственно.
upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
down_part1 = kernel_squared_module * delta_x_pow4;
m_pow2 = m * m;
n_pow2 = n * n;
down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph3.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Две матрицы размера <math>N\times N</math> - спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Матрица размера <math>N\times N</math> - спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.

== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==
Все запуски программы проводились на с/к Ломоносов 2. Для распараллеливания использовалась технология OpenMP.

{| class="wikitable sortable"
|+ Характеристики системы
|-
! Характеристика || Значение
|-
| Пиковая производительность || 5,505 Пфлопс
|-
| Производительность на тесте Linpack || 2,478 Пфлопс
|-
| Число вычислительных узлов || 1 730
|-
| Основной тип процессора || Intel Haswell-EP E5-2697v3, 2.6 GHz, 14 cores
|-
| Тип ускорителя || NVidia Tesla K40M
|-
| Общее число ядер || 64 384
|-
| Оперативная память на узел || 64 GB
|-
| Проектная мощность || 10 МВт
|-
| Операционная система || CentOS 7
|}

Конфигурация системы в разделе compute (в котором производились запуски программы):
{| class="wikitable sortable"
|+ 1487 узлов со следующими характеристиками
|-
! Характеристика || Значение
|-
| Объём памяти на GPU || 11.56 GB
|-
| Центральный процессор || Intel Haswell-EP E5-2697v3, 2.6 GHz, 14 cores
|-
| Графический ускоритель || NVidia Tesla K40M
|-
| ОЗУ || 64 GB
|}

Перейдём к графику времени работы алгоритма. На графике по оси Ox, Oy откладывается длина одной стороны изображения в пикселях (т.е. это - число <math>N</math>, если всё изображение - матрица размера <math>N\times N</math>) и количество нитей OpenMP (нити запускаются параллельно, если доступно необходимое количество ресурсов). По оси Oz откладывается время работы в секундах для каждой точки <math>(x, y)</math>. Временем работы здесь является время работы функции, применяющей метод регуляризации Тихонова к изображению и ядру и заполняющей выходную матрицу. Время чтения/записи файлов, а также время выполнения прочих служебных функций не учитывалось.
Количество нитей задаётся следующей последовательностью: <math>1, 4, 8, 12, 16, 20, ..., 128</math>. А длина стороны изображения в пикселях задаётся такой последовательностью: <math>128, 128+16, ..., 256, 256+32, ..., 4096</math>
[[Файл:TikhonovRegularizationTimes4.png|center|'''Рисунок 2.''' Время выполнения параллельной реализации метода регуляризации Тихонова в секундах.]]

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==
Автору такие не известны. Поиск в открытых источниках результата не дал.

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Участник:Kiseliov/Метод регуляризации Тихонова

2022-11-03T10:28:53Z

Kiseliov: /* Масштабируемость алгоритма и его реализации */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительным ядром являются вычисления по формуле из пункта 1.2 в позиции <math>(i, j)</math>, поскольку формула применяется к позиции <math>(i, j)</math> независимо от других точек.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на псевдокоде в стиле C:

<syntaxhighlight lang="C">
tikhonov_regularization(N, *image, *kernel, R, alpha, r, *restored)
{
// ДельтаX в формуле (логично посчитать эти значения единожды)
delta_x = (4 * R) / N;
delta_x_pow2 = delta_x * delta_x;
delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Это выражение - константа, и его логично посчитать заранее. CONSTANT_PI - число пи
constant = (CONSTANT_PI / (2 * R)) * (CONSTANT_PI / (2 * R));

// Цикл по позициям (i, j)
for (i = 0; i < N; i++)
for (j = 0; j < N; j++)
{
// Предполагается, что матрицы хранятся в виде 1D массива, поэтому нужно вычислить правильный индекс
index = i * N + j;

// Мы сдвигаем центр системы координат в центр изображения
m = i - N / 2;
n = N / 2 - j;

// Вычисление формулы. Предполагается, что pow возводит число в степень, conj возвращает сопряжённое комплексное число,
// .real() и .imag() - получение действительной или мнимой частей комплексного числа соответственно.
upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
down_part1 = kernel_squared_module * delta_x_pow4;
m_pow2 = m * m;
n_pow2 = n * n;
down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph3.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Две матрицы размера <math>N\times N</math> - спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Матрица размера <math>N\times N</math> - спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.

== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==
Все запуски программы проводились на с/к Ломоносов 2. Для распараллеливания использовалась технология OpenMP.

{| class="wikitable sortable"
|+ Характеристики системы
|-
! Характеристика || Значение
|-
| Пиковая производительность || 5,505 Пфлопс
|-
| Производительность на тесте Linpack || 2,478 Пфлопс
|-
| Число вычислительных узлов || 1 730
|-
| Основной тип процессора || Intel Haswell-EP E5-2697v3, 2.6 GHz, 14 cores
|-
| Тип ускорителя || NVidia Tesla K40M
|-
| Общее число ядер || 64 384
|-
| Оперативная память на узел || 64 GB
|-
| Проектная мощность || 10 МВт
|-
| Операционная система || CentOS 7
|}

Конфигурация системы в разделе compute (в котором производились запуски программы):
{| class="wikitable sortable"
|+ 1487 узлов со следующими характеристиками
|-
! Характеристика || Значение
|-
| Объём памяти на GPU || 11.56 GB
|-
| Центральный процессор || Intel Haswell-EP E5-2697v3, 2.6 GHz, 14 cores
|-
| Графический ускоритель || NVidia Tesla K40M
|-
| ОЗУ || 64 GB
|}

Перейдём к графику времени работы алгоритма. На графике по оси Ox, Oy откладывается длина одной стороны изображения в пикселях (т.е. это - число <math>N</math>, если всё изображение - матрица размера <math>N\times N</math>) и количество нитей OpenMP (нити запускаются параллельно, если доступно необходимое количество ресурсов). По оси Oz откладывается время работы в секундах для каждой точки <math>(x, y)</math>. Временем работы здесь является время работы функции, применяющей метод регуляризации Тихонова к изображению и ядру и заполняющей выходную матрицу. Время чтения/записи файлов, а также время выполнения прочих служебных функций не учитывалось.
Количество нитей задаётся следующей последовательностью: <math>1, 4, 8, 12, 16, 20, ..., 128</math>. А длина стороны изображения в пикселях задаётся такой последовательностью: <math>128, 128+16, ..., 256, 256+32, ..., 4096</math>
[[Файл:TikhonovRegularizationTimes4.png|center|'''Рисунок 2.''' Время выполнения параллельной реализации метода регуляризации Тихонова в секундах.]]

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Участник:Kiseliov/Метод регуляризации Тихонова

2022-11-03T10:23:57Z

Kiseliov: /* Масштабируемость алгоритма и его реализации */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительным ядром являются вычисления по формуле из пункта 1.2 в позиции <math>(i, j)</math>, поскольку формула применяется к позиции <math>(i, j)</math> независимо от других точек.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на псевдокоде в стиле C:

<syntaxhighlight lang="C">
tikhonov_regularization(N, *image, *kernel, R, alpha, r, *restored)
{
// ДельтаX в формуле (логично посчитать эти значения единожды)
delta_x = (4 * R) / N;
delta_x_pow2 = delta_x * delta_x;
delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Это выражение - константа, и его логично посчитать заранее. CONSTANT_PI - число пи
constant = (CONSTANT_PI / (2 * R)) * (CONSTANT_PI / (2 * R));

// Цикл по позициям (i, j)
for (i = 0; i < N; i++)
for (j = 0; j < N; j++)
{
// Предполагается, что матрицы хранятся в виде 1D массива, поэтому нужно вычислить правильный индекс
index = i * N + j;

// Мы сдвигаем центр системы координат в центр изображения
m = i - N / 2;
n = N / 2 - j;

// Вычисление формулы. Предполагается, что pow возводит число в степень, conj возвращает сопряжённое комплексное число,
// .real() и .imag() - получение действительной или мнимой частей комплексного числа соответственно.
upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
down_part1 = kernel_squared_module * delta_x_pow4;
m_pow2 = m * m;
n_pow2 = n * n;
down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph3.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Две матрицы размера <math>N\times N</math> - спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Матрица размера <math>N\times N</math> - спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.

== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==
Все запуски программы проводились на с/к Ломоносов 2. Для распараллеливания использовалась технология OpenMP.

{| class="wikitable sortable"
|+ Характеристики системы
|-
! Характеристика || Значение
|-
| Пиковая производительность || 5,505 Пфлопс
|-
| Производительность на тесте Linpack || 2,478 Пфлопс
|-
| Число вычислительных узлов || 1 730
|-
| Основной тип процессора || Intel Haswell-EP E5-2697v3, 2.6 GHz, 14 cores
|-
| Тип ускорителя || NVidia Tesla K40M
|-
| Общее число ядер || 64 384
|-
| Оперативная память на узел || 64 GB
|-
| Проектная мощность || 10 МВт
|-
| Операционная система || CentOS 7
|}

Конфигурация системы в разделе compute (в котором производились запуски программы):
{| class="wikitable sortable"
|+ 1487 узлов со следующими характеристиками
|-
! Характеристика || Значение
|-
| Объём памяти на GPU || 11.56 GB
|-
| Центральный процессор || Intel Haswell-EP E5-2697v3, 2.6 GHz, 14 cores
|-
| Графический ускоритель || NVidia Tesla K40M
|-
| ОЗУ || 64 GB
|}

Перейдём к графику времени работы алгоритма. На графике по оси Ox, Oy откладывается длина одной стороны изображения в пикселях (т.е. это - число <math>N</math>, если всё изображение - матрица размера <math>N\times N</math>) и количество нитей OpenMP (нити запускаются параллельно, если доступно необходимое количество ресурсов). По оси Oz откладывается время работы в секундах для каждой точки <math>(x, y)</math>. Временем работы здесь является время работы функции, применяющей метод регуляризации Тихонова к изображению и ядру и заполняющей выходную матрицу. Время чтения/записи файлов, а также время выполнения прочих служебных функций не учитывалось.
Количество нитей задаётся степенью двойки. Длина стороны изображения в пикселях задаётся следующей последовательностью: <math>128, 128+16, ..., 256, 256+32, ..., 4096</math>
[[Файл:TikhonovRegularizationTimes4.png|center|'''Рисунок 2.''' Время выполнения параллельной реализации метода регуляризации Тихонова в секундах.]]

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Файл:TikhonovRegularizationTimes4.png

2022-11-03T10:21:02Z

Kiseliov:

График времени работы алгоритма

Участник:Kiseliov/Метод регуляризации Тихонова

2022-11-02T14:08:08Z

Kiseliov: /* Масштабируемость алгоритма и его реализации */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительным ядром являются вычисления по формуле из пункта 1.2 в позиции <math>(i, j)</math>, поскольку формула применяется к позиции <math>(i, j)</math> независимо от других точек.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на псевдокоде в стиле C:

<syntaxhighlight lang="C">
tikhonov_regularization(N, *image, *kernel, R, alpha, r, *restored)
{
// ДельтаX в формуле (логично посчитать эти значения единожды)
delta_x = (4 * R) / N;
delta_x_pow2 = delta_x * delta_x;
delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Это выражение - константа, и его логично посчитать заранее. CONSTANT_PI - число пи
constant = (CONSTANT_PI / (2 * R)) * (CONSTANT_PI / (2 * R));

// Цикл по позициям (i, j)
for (i = 0; i < N; i++)
for (j = 0; j < N; j++)
{
// Предполагается, что матрицы хранятся в виде 1D массива, поэтому нужно вычислить правильный индекс
index = i * N + j;

// Мы сдвигаем центр системы координат в центр изображения
m = i - N / 2;
n = N / 2 - j;

// Вычисление формулы. Предполагается, что pow возводит число в степень, conj возвращает сопряжённое комплексное число,
// .real() и .imag() - получение действительной или мнимой частей комплексного числа соответственно.
upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
down_part1 = kernel_squared_module * delta_x_pow4;
m_pow2 = m * m;
n_pow2 = n * n;
down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph3.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Две матрицы размера <math>N\times N</math> - спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Матрица размера <math>N\times N</math> - спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.

== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==
Все запуски программы проводились на с/к Ломоносов 2. Для распараллеливания использовалась технология OpenMP.

{| class="wikitable sortable"
|+ Характеристики системы
|-
! Характеристика || Значение
|-
| Пиковая производительность || 5,505 Пфлопс
|-
| Производительность на тесте Linpack || 2,478 Пфлопс
|-
| Число вычислительных узлов || 1 730
|-
| Основной тип процессора || Intel Haswell-EP E5-2697v3, 2.6 GHz, 14 cores
|-
| Тип ускорителя || NVidia Tesla K40M
|-
| Общее число ядер || 64 384
|-
| Оперативная память на узел || 64 GB
|-
| Проектная мощность || 10 МВт
|-
| Операционная система || CentOS 7
|}

Конфигурация системы в разделе compute (в котором производились запуски программы):
{| class="wikitable sortable"
|+ 1487 узлов со следующими характеристиками
|-
! Характеристика || Значение
|-
| Объём памяти на GPU || 11.56 GB
|-
| Центральный процессор || Intel Haswell-EP E5-2697v3, 2.6 GHz, 14 cores
|-
| Графический ускоритель || NVidia Tesla K40M
|-
| ОЗУ || 64 GB
|}

Перейдём к графику времени работы алгоритма. На графике по оси Ox, Oy откладывается длина одной стороны изображения в пикселях (т.е. это - число <math>N</math>, если всё изображение - матрица размера <math>N\times N</math>) и количество процессоров (именно процессоров; это связано с особенностями технологии OpenMP). По оси Oz откладывается время работы в секундах для каждой точки <math>(x, y)</math>. Временем работы здесь является время работы функции, применяющей метод регуляризации Тихонова к изображению и ядру и заполняющей выходную матрицу. Время чтения/записи файлов, а также время выполнения прочих служебных функций не учитывалось.
Количество процессоров задаётся степенью двойки. Длина стороны изображения в пикселях задаётся следующей последовательностью: <math>128, 128+16, ..., 256, 256+32, ..., 4096</math>
[[Файл:TikhonovRegularizationTimes3.png|center|'''Рисунок 2.''' Время выполнения параллельной реализации метода регуляризации Тихонова в секундах.]]

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Файл:TikhonovRegularizationTimes3.png

2022-11-02T14:07:01Z

Kiseliov:

График времени работы алгоритма

Участник:Kiseliov/Метод регуляризации Тихонова

2022-11-02T14:05:37Z

Kiseliov: /* Масштабируемость алгоритма и его реализации */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительным ядром являются вычисления по формуле из пункта 1.2 в позиции <math>(i, j)</math>, поскольку формула применяется к позиции <math>(i, j)</math> независимо от других точек.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на псевдокоде в стиле C:

<syntaxhighlight lang="C">
tikhonov_regularization(N, *image, *kernel, R, alpha, r, *restored)
{
// ДельтаX в формуле (логично посчитать эти значения единожды)
delta_x = (4 * R) / N;
delta_x_pow2 = delta_x * delta_x;
delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Это выражение - константа, и его логично посчитать заранее. CONSTANT_PI - число пи
constant = (CONSTANT_PI / (2 * R)) * (CONSTANT_PI / (2 * R));

// Цикл по позициям (i, j)
for (i = 0; i < N; i++)
for (j = 0; j < N; j++)
{
// Предполагается, что матрицы хранятся в виде 1D массива, поэтому нужно вычислить правильный индекс
index = i * N + j;

// Мы сдвигаем центр системы координат в центр изображения
m = i - N / 2;
n = N / 2 - j;

// Вычисление формулы. Предполагается, что pow возводит число в степень, conj возвращает сопряжённое комплексное число,
// .real() и .imag() - получение действительной или мнимой частей комплексного числа соответственно.
upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
down_part1 = kernel_squared_module * delta_x_pow4;
m_pow2 = m * m;
n_pow2 = n * n;
down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph3.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Две матрицы размера <math>N\times N</math> - спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Матрица размера <math>N\times N</math> - спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.

== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==
Все запуски программы проводились на с/к Ломоносов 2. Для распараллеливания использовалась технология OpenMP.

{| class="wikitable sortable"
|+ Характеристики системы
|-
! Характеристика || Значение
|-
| Пиковая производительность || 5,505 Пфлопс
|-
| Производительность на тесте Linpack || 2,478 Пфлопс
|-
| Число вычислительных узлов || 1 730
|-
| Основной тип процессора || Intel Haswell-EP E5-2697v3, 2.6 GHz, 14 cores
|-
| Тип ускорителя || NVidia Tesla K40M
|-
| Общее число ядер || 64 384
|-
| Оперативная память на узел || 64 GB
|-
| Проектная мощность || 10 МВт
|-
| Операционная система || CentOS 7
|}

Конфигурация системы в разделе compute (в котором производились запуски программы):
{| class="wikitable sortable"
|+ 1487 узлов со следующими характеристиками
|-
! Характеристика || Значение
|-
| Объём памяти на GPU || 11.56 GB
|-
| Центральный процессор || Intel Haswell-EP E5-2697v3, 2.6 GHz, 14 cores
|-
| Графический ускоритель || NVidia Tesla K40M
|-
| ОЗУ || 64 GB
|}

Перейдём к графику времени работы алгоритма. На графике по оси Ox, Oy откладывается длина одной стороны изображения в пикселях (т.е. это - число <math>N</math>, если всё изображение - матрица размера <math>N\times N</math>) и количество процессоров (именно процессоров; это связано с особенностями технологии OpenMP). По оси Oz откладывается время работы в секундах для каждой точки <math>(x, y)</math>. Временем работы здесь является время работы функции, применяющей метод регуляризации Тихонова к изображению и ядру и заполняющей выходную матрицу. Время чтения/записи файлов, а также время выполнения прочих служебных функций не учитывалось.
Количество процессоров задаётся степенью двойки. Длина стороны изображения в пикселях задаётся следующей последовательностью: <math>128, 128+16, ..., 256, 256+32, ..., 4096</math>
[[Файл:TikhonovRegularizationTimes2.png|center|'''Рисунок 2.''' Время выполнения параллельной реализации метода регуляризации Тихонова в секундах.]]

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Участник:Kiseliov/Метод регуляризации Тихонова

2022-11-02T14:04:02Z

Kiseliov: /* Масштабируемость алгоритма и его реализации */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительным ядром являются вычисления по формуле из пункта 1.2 в позиции <math>(i, j)</math>, поскольку формула применяется к позиции <math>(i, j)</math> независимо от других точек.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на псевдокоде в стиле C:

<syntaxhighlight lang="C">
tikhonov_regularization(N, *image, *kernel, R, alpha, r, *restored)
{
// ДельтаX в формуле (логично посчитать эти значения единожды)
delta_x = (4 * R) / N;
delta_x_pow2 = delta_x * delta_x;
delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Это выражение - константа, и его логично посчитать заранее. CONSTANT_PI - число пи
constant = (CONSTANT_PI / (2 * R)) * (CONSTANT_PI / (2 * R));

// Цикл по позициям (i, j)
for (i = 0; i < N; i++)
for (j = 0; j < N; j++)
{
// Предполагается, что матрицы хранятся в виде 1D массива, поэтому нужно вычислить правильный индекс
index = i * N + j;

// Мы сдвигаем центр системы координат в центр изображения
m = i - N / 2;
n = N / 2 - j;

// Вычисление формулы. Предполагается, что pow возводит число в степень, conj возвращает сопряжённое комплексное число,
// .real() и .imag() - получение действительной или мнимой частей комплексного числа соответственно.
upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
down_part1 = kernel_squared_module * delta_x_pow4;
m_pow2 = m * m;
n_pow2 = n * n;
down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph3.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Две матрицы размера <math>N\times N</math> - спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Матрица размера <math>N\times N</math> - спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.

== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==
Все запуски программы проводились на с/к Ломоносов 2. Для распараллеливания использовалась технология OpenMP.

{| class="wikitable sortable"
|+ Характеристики системы
|-
! Характеристика || Значение
|-
| Пиковая производительность || 5,505 Пфлопс
|-
| Производительность на тесте Linpack || 2,478 Пфлопс
|-
| Число вычислительных узлов || 1 730
|-
| Основной тип процессора || Intel Haswell-EP E5-2697v3, 2.6 GHz, 14 cores
|-
| Тип ускорителя || NVidia Tesla K40M
|-
| Общее число ядер || 64 384
|-
| Оперативная память на узел || 64 GB
|-
| Проектная мощность || 10 МВт
|-
| Операционная система || CentOS 7
|}

Конфигурация системы в разделе compute (в котором производились запуски программы):
{| class="wikitable sortable"
|+ 1487 узлов со следующими характеристиками
|-
! Характеристика || Значение
|-
| Объём памяти на GPU || 11.56 GB
|-
| Центральный процессор || Intel Haswell-EP E5-2697v3, 2.6 GHz, 14 cores
|-
| Графический ускоритель || NVidia Tesla K40M
|-
| ОЗУ || 64 GB
|}

Перейдём к графику времени работы алгоритма. На графике по оси Ox, Oy откладывается размер одной стороны изображения в пикселях (т.е. это - число <math>N</math>, если всё изображение - матрица размера <math>N\times N</math>) и количество процессоров (именно процессоров; Это связано с особенностями технологии OpenMP). По оси Oz откладывается время работы в секундах для каждой точки <math>(x, y)</math>. Временем работы здесь является время работы функции, выполняющей метод регуляризации Тихонова над изображением и ядром и заполняющей выходную матрицу. Время чтения/записи файлов, а также время выполнения прочих служебных функций не учитывалось.
Количество процессоров задаётся степенью двойки. Длина стороны изображения в пикселях задаётся следующей последовательностью: <math>128, 128+16, ..., 256, 256+32, ..., 4096</math>
[[Файл:TikhonovRegularizationTimes2.png|center|'''Рисунок 2.''' Время выполнения параллельной реализации метода регуляризации Тихонова в секундах.]]

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Файл:TikhonovRegularizationTimes2.png

2022-11-02T14:03:40Z

Kiseliov:

График времени работы алгоритма

Участник:Kiseliov/Метод регуляризации Тихонова

2022-11-02T13:14:28Z

Kiseliov: /* Масштабируемость алгоритма и его реализации */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительным ядром являются вычисления по формуле из пункта 1.2 в позиции <math>(i, j)</math>, поскольку формула применяется к позиции <math>(i, j)</math> независимо от других точек.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на псевдокоде в стиле C:

<syntaxhighlight lang="C">
tikhonov_regularization(N, *image, *kernel, R, alpha, r, *restored)
{
// ДельтаX в формуле (логично посчитать эти значения единожды)
delta_x = (4 * R) / N;
delta_x_pow2 = delta_x * delta_x;
delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Это выражение - константа, и его логично посчитать заранее. CONSTANT_PI - число пи
constant = (CONSTANT_PI / (2 * R)) * (CONSTANT_PI / (2 * R));

// Цикл по позициям (i, j)
for (i = 0; i < N; i++)
for (j = 0; j < N; j++)
{
// Предполагается, что матрицы хранятся в виде 1D массива, поэтому нужно вычислить правильный индекс
index = i * N + j;

// Мы сдвигаем центр системы координат в центр изображения
m = i - N / 2;
n = N / 2 - j;

// Вычисление формулы. Предполагается, что pow возводит число в степень, conj возвращает сопряжённое комплексное число,
// .real() и .imag() - получение действительной или мнимой частей комплексного числа соответственно.
upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
down_part1 = kernel_squared_module * delta_x_pow4;
m_pow2 = m * m;
n_pow2 = n * n;
down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph3.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Две матрицы размера <math>N\times N</math> - спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Матрица размера <math>N\times N</math> - спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.

== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==
Все запуски программы проводились на с/к Ломоносов 2. Для распараллеливания использовалась технология OpenMP.

{| class="wikitable sortable"
|+ Характеристики системы
|-
! Характеристика || Значение
|-
| Пиковая производительность || 5,505 Пфлопс
|-
| Производительность на тесте Linpack || 2,478 Пфлопс
|-
| Число вычислительных узлов || 1 730
|-
| Основной тип процессора || Intel Haswell-EP E5-2697v3, 2.6 GHz, 14 cores
|-
| Тип ускорителя || NVidia Tesla K40M
|-
| Общее число ядер || 64 384
|-
| Оперативная память на узел || 64 GB
|-
| Проектная мощность || 10 МВт
|-
| Операционная система || CentOS 7
|}

Конфигурация системы в разделе compute (в котором производились запуски программы):
{| class="wikitable sortable"
|+ 1487 узлов со следующими характеристиками
|-
! Характеристика || Значение
|-
| Объём памяти на GPU || 11.56 GB
|-
| Центральный процессор || Intel Haswell-EP E5-2697v3, 2.6 GHz, 14 cores
|-
| Графический ускоритель || NVidia Tesla K40M
|-
| ОЗУ || 64 GB
|}

Перейдём к замерам времени работы алгоритма.
[[Файл:TikhonovRegularizationTimes1.png|center|'''Рисунок 2.''' Время выполнения параллельной реализации метода регуляризации Тихонова в секундах.]]

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-27T18:22:24Z

Kiseliov: /* Масштабируемость алгоритма и его реализации */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительным ядром являются вычисления по формуле из пункта 1.2 в позиции <math>(i, j)</math>, поскольку формула применяется к позиции <math>(i, j)</math> независимо от других точек.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на псевдокоде в стиле C:

<syntaxhighlight lang="C">
tikhonov_regularization(N, *image, *kernel, R, alpha, r, *restored)
{
// ДельтаX в формуле (логично посчитать эти значения единожды)
delta_x = (4 * R) / N;
delta_x_pow2 = delta_x * delta_x;
delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Это выражение - константа, и его логично посчитать заранее. CONSTANT_PI - число пи
constant = (CONSTANT_PI / (2 * R)) * (CONSTANT_PI / (2 * R));

// Цикл по позициям (i, j)
for (i = 0; i < N; i++)
for (j = 0; j < N; j++)
{
// Предполагается, что матрицы хранятся в виде 1D массива, поэтому нужно вычислить правильный индекс
index = i * N + j;

// Мы сдвигаем центр системы координат в центр изображения
m = i - N / 2;
n = N / 2 - j;

// Вычисление формулы. Предполагается, что pow возводит число в степень, conj возвращает сопряжённое комплексное число,
// .real() и .imag() - получение действительной или мнимой частей комплексного числа соответственно.
upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
down_part1 = kernel_squared_module * delta_x_pow4;
m_pow2 = m * m;
n_pow2 = n * n;
down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph3.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Две матрицы размера <math>N\times N</math> - спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Матрица размера <math>N\times N</math> - спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.

== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==
Все запуски программы проводились на с/к Ломоносов 2. Для распараллеливания использовалась технология OpenMP.

{| class="wikitable sortable"
|+ Характеристики системы
|-
! Характеристика || Значение
|-
| Пиковая производительность || 2,962 Пфлопс
|-
| Производительность на тесте Linpack || 2,1 Пфлопс
|-
| Число вычислительных узлов || 1 472
|-
| Основной тип процессора || Intel Haswell-EP E5-2697v3, 2.6 GHz, 14 cores
|-
| Тип ускорителя || NVidia Tesla K40M
|-
| Общее число ядер || 42 688
|-
| Оперативная память на узел || 64 GB
|-
| Оперативная память на стойку || 16 Тбайт
|-
| Проектная мощность || 10 МВт
|-
| Расчетная производительность || 2,5 Пфлопс
|-
| Максимальное число процессоров в стойке || 256
|-
| Максимальное число графических ускорителей в стойке || 256
|-
| Операционная система || CentOS 7
|}

Конфигурация системы в разделе compute (в котором производились запуски программы):
{| class="wikitable sortable"
|+ 1504 узлов со следующими характеристиками
|-
! Характеристика || Значение
|-
| Объём памяти на GPU || 11.56 GB
|-
| GPU || 1x Tesla K40s
|-
| CPU || 1x Intel Xeon E5-2697 v3 2.60GHz
|-
| CPU-ядер || 14
|-
| ОЗУ || 64 GB
|}

Перейдём к замерам времени работы алгоритма.
[[Файл:TikhonovRegularizationTimes1.png|center|'''Рисунок 2.''' Время выполнения параллельной реализации метода регуляризации Тихонова в секундах.]]

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Файл:TikhonovRegularizationTimes1.png

2022-10-27T18:19:33Z

Kiseliov:

График времени работы алгоритма

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-26T19:47:56Z

Kiseliov: /* Масштабируемость алгоритма и его реализации */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительным ядром являются вычисления по формуле из пункта 1.2 в позиции <math>(i, j)</math>, поскольку формула применяется к позиции <math>(i, j)</math> независимо от других точек.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на псевдокоде в стиле C:

<syntaxhighlight lang="C">
tikhonov_regularization(N, *image, *kernel, R, alpha, r, *restored)
{
// ДельтаX в формуле (логично посчитать эти значения единожды)
delta_x = (4 * R) / N;
delta_x_pow2 = delta_x * delta_x;
delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Это выражение - константа, и его логично посчитать заранее. CONSTANT_PI - число пи
constant = (CONSTANT_PI / (2 * R)) * (CONSTANT_PI / (2 * R));

// Цикл по позициям (i, j)
for (i = 0; i < N; i++)
for (j = 0; j < N; j++)
{
// Предполагается, что матрицы хранятся в виде 1D массива, поэтому нужно вычислить правильный индекс
index = i * N + j;

// Мы сдвигаем центр системы координат в центр изображения
m = i - N / 2;
n = N / 2 - j;

// Вычисление формулы. Предполагается, что pow возводит число в степень, conj возвращает сопряжённое комплексное число,
// .real() и .imag() - получение действительной или мнимой частей комплексного числа соответственно.
upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
down_part1 = kernel_squared_module * delta_x_pow4;
m_pow2 = m * m;
n_pow2 = n * n;
down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph3.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Две матрицы размера <math>N\times N</math> - спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Матрица размера <math>N\times N</math> - спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.

== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==
Все запуски программы проводились на с/к Ломоносов 2. Для распараллеливания использовалась технология OpenMP.

{| class="wikitable sortable"
|+ Характеристики системы
|-
! Характеристика || Значение
|-
| Пиковая производительность || 2,962 Пфлопс
|-
| Производительность на тесте Linpack || 2,1 Пфлопс
|-
| Число вычислительных узлов || 1 472
|-
| Основной тип процессора || Intel Haswell-EP E5-2697v3, 2.6 GHz, 14 cores
|-
| Тип ускорителя || NVidia Tesla K40M
|-
| Общее число ядер || 42 688
|-
| Оперативная память на узел || 64 GB
|-
| Оперативная память на стойку || 16 Тбайт
|-
| Проектная мощность || 10 МВт
|-
| Расчетная производительность || 2,5 Пфлопс
|-
| Максимальное число процессоров в стойке || 256
|-
| Максимальное число графических ускорителей в стойке || 256
|-
| Операционная система || CentOS 7
|}

Конфигурация системы в разделе compute (в котором производились запуски программы):
{| class="wikitable sortable"
|+ 1504 узлов со следующими характеристиками
|-
! Характеристика || Значение
|-
| Объём памяти на GPU || 11.56 GB
|-
| GPU || 1x Tesla K40s
|-
| CPU || 1x Intel Xeon E5-2697 v3 2.60GHz
|-
| CPU-ядер || 14
|-
| ОЗУ || 64 GB
|}

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-25T18:03:55Z

Kiseliov: /* Схема реализации последовательного алгоритма */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительным ядром являются вычисления по формуле из пункта 1.2 в позиции <math>(i, j)</math>, поскольку формула применяется к позиции <math>(i, j)</math> независимо от других точек.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на псевдокоде в стиле C:

<syntaxhighlight lang="C">
tikhonov_regularization(N, *image, *kernel, R, alpha, r, *restored)
{
// ДельтаX в формуле (логично посчитать эти значения единожды)
delta_x = (4 * R) / N;
delta_x_pow2 = delta_x * delta_x;
delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Это выражение - константа, и его логично посчитать заранее. CONSTANT_PI - число пи
constant = (CONSTANT_PI / (2 * R)) * (CONSTANT_PI / (2 * R));

// Цикл по позициям (i, j)
for (i = 0; i < N; i++)
for (j = 0; j < N; j++)
{
// Предполагается, что матрицы хранятся в виде 1D массива, поэтому нужно вычислить правильный индекс
index = i * N + j;

// Мы сдвигаем центр системы координат в центр изображения
m = i - N / 2;
n = N / 2 - j;

// Вычисление формулы. Предполагается, что pow возводит число в степень, conj возвращает сопряжённое комплексное число,
// .real() и .imag() - получение действительной или мнимой частей комплексного числа соответственно.
upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
down_part1 = kernel_squared_module * delta_x_pow4;
m_pow2 = m * m;
n_pow2 = n * n;
down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph3.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Две матрицы размера <math>N\times N</math> - спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Матрица размера <math>N\times N</math> - спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.

== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-21T17:14:43Z

Kiseliov: /* Входные и выходные данные алгоритма */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительным ядром являются вычисления по формуле из пункта 1.2 в позиции <math>(i, j)</math>, поскольку формула применяется к позиции <math>(i, j)</math> независимо от других точек.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на C++:

<syntaxhighlight lang="cpp">
#define _USE_MATH_DEFINES

#include <cmath>
#include <complex>
#include <complex.h>

using namespace std;

int TikhonovReg(const int N, const complex<long double>* image, const complex<long double>* kernel, const int R, const double alpha, const int r, complex<long double>* restored)
{
// DeltaX in formula (we want to compute this values once)
double delta_x = (4 * R) / double(N);
long double delta_x_pow2 = delta_x * delta_x;
long double delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Constant part in formula, that can be precomputed
long double constant = (M_PI / (2 * R)) * (M_PI / (2 * R));

// For loop gets parallel. Each step of loop for point (i, j) is executed separately
for (int i = 0; i < N; i++)
for (int j = 0; j < N; j++)
{
// Index (i, j) in 1d mass
long index = i * N + j;

// Move coordonate system center to the center of image
int m = i - N / 2;
int n = N / 2 - j;

// Calculating formula
complex<long double> upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
complex<long double> kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
complex<long double> down_part1 = kernel_squared_module * delta_x_pow4;
int m_pow2 = m * m;
int n_pow2 = n * n;
complex<long double> down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
complex<long double>down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph3.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Две матрицы размера <math>N\times N</math> - спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Матрица размера <math>N\times N</math> - спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.

== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-21T17:13:52Z

Kiseliov: /* Входные и выходные данные алгоритма */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительным ядром являются вычисления по формуле из пункта 1.2 в позиции <math>(i, j)</math>, поскольку формула применяется к позиции <math>(i, j)</math> независимо от других точек.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на C++:

<syntaxhighlight lang="cpp">
#define _USE_MATH_DEFINES

#include <cmath>
#include <complex>
#include <complex.h>

using namespace std;

int TikhonovReg(const int N, const complex<long double>* image, const complex<long double>* kernel, const int R, const double alpha, const int r, complex<long double>* restored)
{
// DeltaX in formula (we want to compute this values once)
double delta_x = (4 * R) / double(N);
long double delta_x_pow2 = delta_x * delta_x;
long double delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Constant part in formula, that can be precomputed
long double constant = (M_PI / (2 * R)) * (M_PI / (2 * R));

// For loop gets parallel. Each step of loop for point (i, j) is executed separately
for (int i = 0; i < N; i++)
for (int j = 0; j < N; j++)
{
// Index (i, j) in 1d mass
long index = i * N + j;

// Move coordonate system center to the center of image
int m = i - N / 2;
int n = N / 2 - j;

// Calculating formula
complex<long double> upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
complex<long double> kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
complex<long double> down_part1 = kernel_squared_module * delta_x_pow4;
int m_pow2 = m * m;
int n_pow2 = n * n;
complex<long double> down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
complex<long double>down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph3.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Две матрицы размера - спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Матрица размера - спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.

== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-19T12:05:37Z

Kiseliov: /* Информационный граф */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительным ядром являются вычисления по формуле из пункта 1.2 в позиции <math>(i, j)</math>, поскольку формула применяется к позиции <math>(i, j)</math> независимо от других точек.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на C++:

<syntaxhighlight lang="cpp">
#define _USE_MATH_DEFINES

#include <cmath>
#include <complex>
#include <complex.h>

using namespace std;

int TikhonovReg(const int N, const complex<long double>* image, const complex<long double>* kernel, const int R, const double alpha, const int r, complex<long double>* restored)
{
// DeltaX in formula (we want to compute this values once)
double delta_x = (4 * R) / double(N);
long double delta_x_pow2 = delta_x * delta_x;
long double delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Constant part in formula, that can be precomputed
long double constant = (M_PI / (2 * R)) * (M_PI / (2 * R));

// For loop gets parallel. Each step of loop for point (i, j) is executed separately
for (int i = 0; i < N; i++)
for (int j = 0; j < N; j++)
{
// Index (i, j) in 1d mass
long index = i * N + j;

// Move coordonate system center to the center of image
int m = i - N / 2;
int n = N / 2 - j;

// Calculating formula
complex<long double> upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
complex<long double> kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
complex<long double> down_part1 = kernel_squared_module * delta_x_pow4;
int m_pow2 = m * m;
int n_pow2 = n * n;
complex<long double> down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
complex<long double>down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph3.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.
== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Файл:TikhonovRegularizationGraph3.png

2022-10-19T12:05:15Z

Kiseliov:

Граф алгоритма

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-19T11:25:33Z

Kiseliov: /* Информационный граф */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительным ядром являются вычисления по формуле из пункта 1.2 в позиции <math>(i, j)</math>, поскольку формула применяется к позиции <math>(i, j)</math> независимо от других точек.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на C++:

<syntaxhighlight lang="cpp">
#define _USE_MATH_DEFINES

#include <cmath>
#include <complex>
#include <complex.h>

using namespace std;

int TikhonovReg(const int N, const complex<long double>* image, const complex<long double>* kernel, const int R, const double alpha, const int r, complex<long double>* restored)
{
// DeltaX in formula (we want to compute this values once)
double delta_x = (4 * R) / double(N);
long double delta_x_pow2 = delta_x * delta_x;
long double delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Constant part in formula, that can be precomputed
long double constant = (M_PI / (2 * R)) * (M_PI / (2 * R));

// For loop gets parallel. Each step of loop for point (i, j) is executed separately
for (int i = 0; i < N; i++)
for (int j = 0; j < N; j++)
{
// Index (i, j) in 1d mass
long index = i * N + j;

// Move coordonate system center to the center of image
int m = i - N / 2;
int n = N / 2 - j;

// Calculating formula
complex<long double> upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
complex<long double> kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
complex<long double> down_part1 = kernel_squared_module * delta_x_pow4;
int m_pow2 = m * m;
int n_pow2 = n * n;
complex<long double> down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
complex<long double>down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph2.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.
== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Файл:TikhonovRegularizationGraph2.png

2022-10-19T11:25:07Z

Kiseliov:

Граф алгоритма

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-19T11:12:44Z

Kiseliov: /* Схема реализации последовательного алгоритма */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительным ядром являются вычисления по формуле из пункта 1.2 в позиции <math>(i, j)</math>, поскольку формула применяется к позиции <math>(i, j)</math> независимо от других точек.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на C++:

<syntaxhighlight lang="cpp">
#define _USE_MATH_DEFINES

#include <cmath>
#include <complex>
#include <complex.h>

using namespace std;

int TikhonovReg(const int N, const complex<long double>* image, const complex<long double>* kernel, const int R, const double alpha, const int r, complex<long double>* restored)
{
// DeltaX in formula (we want to compute this values once)
double delta_x = (4 * R) / double(N);
long double delta_x_pow2 = delta_x * delta_x;
long double delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Constant part in formula, that can be precomputed
long double constant = (M_PI / (2 * R)) * (M_PI / (2 * R));

// For loop gets parallel. Each step of loop for point (i, j) is executed separately
for (int i = 0; i < N; i++)
for (int j = 0; j < N; j++)
{
// Index (i, j) in 1d mass
long index = i * N + j;

// Move coordonate system center to the center of image
int m = i - N / 2;
int n = N / 2 - j;

// Calculating formula
complex<long double> upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
complex<long double> kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
complex<long double> down_part1 = kernel_squared_module * delta_x_pow4;
int m_pow2 = m * m;
int n_pow2 = n * n;
complex<long double> down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
complex<long double>down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.
== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-19T11:12:26Z

Kiseliov: /* Схема реализации последовательного алгоритма */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительным ядром являются вычисления по формуле из пункта 1.2 в позиции <math>(i, j)</math>, поскольку формула применяется к позиции <math>(i, j)</math> независимо от других точек.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на C++:

<syntaxhighlight lang="cpp">
#define _USE_MATH_DEFINES

#include <cmath>
#include <complex>
#include <complex.h>
#include <omp.h>

using namespace std;

int TikhonovReg(const int N, const complex<long double>* image, const complex<long double>* kernel, const int R, const double alpha, const int r, complex<long double>* restored)
{
// DeltaX in formula (we want to compute this values once)
double delta_x = (4 * R) / double(N);
long double delta_x_pow2 = delta_x * delta_x;
long double delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Constant part in formula, that can be precomputed
long double constant = (M_PI / (2 * R)) * (M_PI / (2 * R));

// For loop gets parallel. Each step of loop for point (i, j) is executed separately
for (int i = 0; i < N; i++)
for (int j = 0; j < N; j++)
{
// Index (i, j) in 1d mass
long index = i * N + j;

// Move coordonate system center to the center of image
int m = i - N / 2;
int n = N / 2 - j;

// Calculating formula
complex<long double> upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
complex<long double> kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
complex<long double> down_part1 = kernel_squared_module * delta_x_pow4;
int m_pow2 = m * m;
int n_pow2 = n * n;
complex<long double> down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
complex<long double>down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.
== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-19T10:49:39Z

Kiseliov: /* Вычислительное ядро алгоритма */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительным ядром являются вычисления по формуле из пункта 1.2 в позиции <math>(i, j)</math>, поскольку формула применяется к позиции <math>(i, j)</math> независимо от других точек.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на C++:

<syntaxhighlight lang="cpp">
#define _USE_MATH_DEFINES

#include <cmath>
#include <complex>
#include <complex.h>
#include <omp.h>

using namespace std;

int TikhonovReg(const int N, const complex<long double>* image, const complex<long double>* kernel, const int R, const double alpha, const int r, complex<long double>* restored)
{
// DeltaX in formula (we want to compute this values once)
double delta_x = (4 * R) / double(N);
long double delta_x_pow2 = delta_x * delta_x;
long double delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Constant part in formula, that can be precomputed
long double constant = (M_PI / (2 * R)) * (M_PI / (2 * R));

// For loop gets parallel. Each step of loop for point (i, j) is executed separately
#pragma omp parallel for collapse(2) num_threads(8)
for (int i = 0; i < N; i++)
for (int j = 0; j < N; j++)
{
// Index (i, j) in 1d mass
long index = i * N + j;

// Move coordonate system center to the center of image
int m = i - N / 2;
int n = N / 2 - j;

// Calculating formula
complex<long double> upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
complex<long double> kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
complex<long double> down_part1 = kernel_squared_module * delta_x_pow4;
int m_pow2 = m * m;
int n_pow2 = n * n;
complex<long double> down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
complex<long double>down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.
== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-18T11:38:26Z

Kiseliov: /* Информационный граф */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительное ядро совпадает с алгоритмом, так как в позиции <math>(i, j)</math> вычисления производятся независимо по формуле, указанной пункте 1.2.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на C++:

<syntaxhighlight lang="cpp">
#define _USE_MATH_DEFINES

#include <cmath>
#include <complex>
#include <complex.h>
#include <omp.h>

using namespace std;

int TikhonovReg(const int N, const complex<long double>* image, const complex<long double>* kernel, const int R, const double alpha, const int r, complex<long double>* restored)
{
// DeltaX in formula (we want to compute this values once)
double delta_x = (4 * R) / double(N);
long double delta_x_pow2 = delta_x * delta_x;
long double delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Constant part in formula, that can be precomputed
long double constant = (M_PI / (2 * R)) * (M_PI / (2 * R));

// For loop gets parallel. Each step of loop for point (i, j) is executed separately
#pragma omp parallel for collapse(2) num_threads(8)
for (int i = 0; i < N; i++)
for (int j = 0; j < N; j++)
{
// Index (i, j) in 1d mass
long index = i * N + j;

// Move coordonate system center to the center of image
int m = i - N / 2;
int n = N / 2 - j;

// Calculating formula
complex<long double> upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
complex<long double> kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
complex<long double> down_part1 = kernel_squared_module * delta_x_pow4;
int m_pow2 = m * m;
int n_pow2 = n * n;
complex<long double> down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
complex<long double>down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода. Рёбра в графе отсутствуют в связи с отсутствием связей по данным.]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.
== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-18T11:34:39Z

Kiseliov: /* Информационный граф */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительное ядро совпадает с алгоритмом, так как в позиции <math>(i, j)</math> вычисления производятся независимо по формуле, указанной пункте 1.2.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на C++:

<syntaxhighlight lang="cpp">
#define _USE_MATH_DEFINES

#include <cmath>
#include <complex>
#include <complex.h>
#include <omp.h>

using namespace std;

int TikhonovReg(const int N, const complex<long double>* image, const complex<long double>* kernel, const int R, const double alpha, const int r, complex<long double>* restored)
{
// DeltaX in formula (we want to compute this values once)
double delta_x = (4 * R) / double(N);
long double delta_x_pow2 = delta_x * delta_x;
long double delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Constant part in formula, that can be precomputed
long double constant = (M_PI / (2 * R)) * (M_PI / (2 * R));

// For loop gets parallel. Each step of loop for point (i, j) is executed separately
#pragma omp parallel for collapse(2) num_threads(8)
for (int i = 0; i < N; i++)
for (int j = 0; j < N; j++)
{
// Index (i, j) in 1d mass
long index = i * N + j;

// Move coordonate system center to the center of image
int m = i - N / 2;
int n = N / 2 - j;

// Calculating formula
complex<long double> upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
complex<long double> kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
complex<long double> down_part1 = kernel_squared_module * delta_x_pow4;
int m_pow2 = m * m;
int n_pow2 = n * n;
complex<long double> down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
complex<long double>down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
[[Файл:TikhonovRegularizationGraph.png|center|thumb|'''Рисунок 1.''' Информационный граф метода регуляризации Тихонова. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)</math> позиции выхода]]

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.
== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-18T11:23:18Z

Kiseliov: /* Информационный граф */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительное ядро совпадает с алгоритмом, так как в позиции <math>(i, j)</math> вычисления производятся независимо по формуле, указанной пункте 1.2.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на C++:

<syntaxhighlight lang="cpp">
#define _USE_MATH_DEFINES

#include <cmath>
#include <complex>
#include <complex.h>
#include <omp.h>

using namespace std;

int TikhonovReg(const int N, const complex<long double>* image, const complex<long double>* kernel, const int R, const double alpha, const int r, complex<long double>* restored)
{
// DeltaX in formula (we want to compute this values once)
double delta_x = (4 * R) / double(N);
long double delta_x_pow2 = delta_x * delta_x;
long double delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Constant part in formula, that can be precomputed
long double constant = (M_PI / (2 * R)) * (M_PI / (2 * R));

// For loop gets parallel. Each step of loop for point (i, j) is executed separately
#pragma omp parallel for collapse(2) num_threads(8)
for (int i = 0; i < N; i++)
for (int j = 0; j < N; j++)
{
// Index (i, j) in 1d mass
long index = i * N + j;

// Move coordonate system center to the center of image
int m = i - N / 2;
int n = N / 2 - j;

// Calculating formula
complex<long double> upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
complex<long double> kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
complex<long double> down_part1 = kernel_squared_module * delta_x_pow4;
int m_pow2 = m * m;
int n_pow2 = n * n;
complex<long double> down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
complex<long double>down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

== Информационный граф ==
В информационном графе отсутствуют связи, так как нет зависимости по данным.
<gallery>
[[Файл:TikhonovRegularizationGraph.png|обрамить|центр|Граф алгоритма. Вершины означают применение формулы из пункта 2 к <math>(i, j)</math> позиции входа для получения <math>(i, j)<math/> позиции выхода)]]
</gallery>

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.
== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Файл:TikhonovRegularizationGraph.png

2022-10-18T11:20:32Z

Kiseliov:

Граф алгоритма

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-11T20:50:28Z

Kiseliov: /* Схема реализации последовательного алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-11T20:48:42Z

Kiseliov: /* Схема реализации последовательного алгоритма */

Автор: Киселёв Е. И.

= ЧАСТЬ. Свойства и структура алгоритмов =
== Общее описание алгоритма ==
Метод регуляризации Тихонова заключается в следующем:

Нам дано некоторое искажённое изображение (в нашем случае мы рассматриваем статические аберрации). Фактически, некоторое изображение искажается при помощи свёртки с так называемым ядром. То есть, мы имеем уравнение Фредгольма первого рода типа свертки вида:

<math>K \circledast z =\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}K(x_1-s_1,x_2-s_2 )z(s_1,s_2 )ds_1\,ds_2 = u(x_1,x_2 ), -\infty<x_1,x_2<\infty</math>

Здесь <math>K(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – аппаратная функция прибора (ядро), <math>u(x_1,x_2 )∈L_2 (\mathbb{R}^2)</math> – искаженное изображение, а <math>z(x_1,x_2 )</math> – искомое реконструируемое изображение.

Наша задача - восстановить исходное изображение, зная параметры ядра. Метод регуляризации Тихонова говорит о том, что решение имеет вид:

<math>z(x_1 ,x_2 ) = \frac{1}{4\pi^2}\textstyle\int\limits_{-\infty}^{\infty}\textstyle\int\limits_{-\infty}^{\infty}\frac{\bar K (-\omega_1,-\omega_2)\bar u (\omega_1,\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}e^{i(\omega_1x_1 + \omega_2x_2)}d\omega_1\,d\omega_2</math>

Здесь <math>\bar K (x_1,x_2 )</math> – спектр ядра, <math>\bar u (x_1,x_2 )</math> – спектр искаженного изображения, а <math>M(\omega_1,\omega_2)</math> – заданная четная функция, обладающая следующими свойствами:
# <math>M(\omega_1,\omega_2)</math> кусочно-непрерывна в любой конечной области
# <math>M(\omega_1,\omega_2)</math> неотрицательна: <math>M(0,0)\ge0</math> и <math>M(\omega_1,\omega_2)>0</math> при <math>\omega_1,\omega_2\neq0</math>
# для достаточно больших <math>|\omega_1|,|\omega_2| \Rightarrow M(\omega_1,\omega_2)\ge C>0</math>
# для <math>\forall \alpha >0 \Rightarrow \frac{\bar K (-\omega_1,-\omega_2)}{|\bar K (\omega_1,\omega_2 )|^2+\alpha M(\omega_1,\omega_2)}∈L_2 (\mathbb{R}^2)</math>

== Математическое описание алгоритма ==
Пусть нам даны спектры искажённого изображения <math>\bar u_{ij}</math> и ядра <math>\bar K_{ij}</math> - это квадратные матрицы, состоящие из комплексных чисел, размера <math>N\times N</math>. Тогда спектр исходного изображения в позиции <math>(i, j)</math> находится по формуле:

<math>\bar z_{ij} =\frac{(\bar K_{ij})^*\bar u_{ij} (\vartriangle x)^2}{\bigl|\bar K_{ij}\bigl|^2(\vartriangle x)^4+\alpha \bigl(\bigl(\frac{\pi}{2R}\bigl)^2(i^2+j^2)\bigl)^r}</math> 

Здесь <math>\vartriangle x=\frac{4R}{N}</math>, где R - один из параметров ядра (известная, наперёд заданная величина).<math>(\bar K_{ij})^*</math> - сопряжение спектра ядра. <math>\alpha</math> и <math>r</math> - параметры метода.

В реализации будем считать <math>\alpha = 0.5</math> и <math>r = 1</math> (наиболее подходящие параметры в нашем случае).

== Вычислительное ядро алгоритма ==
В данном методе вычислительное ядро совпадает с алгоритмом, так как в позиции <math>(i, j)</math> вычисления производятся независимо по формуле, указанной пункте 1.2.

Вычислительная сложность формулы - 3 умножения в числителе, 9 умножений, 2 сложения в знаменателе и одно общее деление (возведение в степень не считаем, так как положили <math>r = 1</math>). Итого, 15 операций. Формулы применяется независимо в каждой точке <math>(i, j)</math> <math>\Rightarrow</math> вычислительная сложность ядра алгоритма - 15 операций.

== Макроструктура алгоритма ==
Макроструктура алгоритма представляет собой проход по всем точкам <math/>(i, j)<math/> и вычисление <math/>(i, j)<math/> позиции результирующей матрицы по формуле, указанной в пункте 1.2.

== Схема реализации последовательного алгоритма ==
Реализация указанного алгоритма на C++:

<syntaxhighlight lang="cpp">
#define _USE_MATH_DEFINES

#include <cmath>
#include <complex>
#include <complex.h>
#include <omp.h>

using namespace std;

int TikhonovReg(const int threads, const int N, const complex<long double>* image, const complex<long double>* kernel, const int R, const double alpha, const int r, complex<long double>* restored)
{
// Set number of threads
omp_set_num_threads(threads);

// DeltaX in formula (we want to compute this values once)
double delta_x = (4 * R) / double(N);
long double delta_x_pow2 = delta_x * delta_x;
long double delta_x_pow4 = delta_x_pow2 * delta_x_pow2;

// Constant part in formula, that can be precomputed
long double constant = (M_PI / (2 * R)) * (M_PI / (2 * R));

// For loop gets parallel. Each step of loop for point (i, j) is executed separately
#pragma omp parallel for collapse(2) num_threads(8)
for (int i = 0; i < N; i++)
for (int j = 0; j < N; j++)
{
// Index (i, j) in 1d mass
long index = i * N + j;

// Move coordonate system center to the center of image
int m = i - N / 2;
int n = N / 2 - j;

// Calculating formula
complex<long double> upper_part = conj(kernel[index]) * image[index] * delta_x_pow2;
complex<long double> kernel_squared_module = kernel[index].real() * kernel[index].real() + kernel[index].imag() * kernel[index].imag();
complex<long double> down_part1 = kernel_squared_module * delta_x_pow4;
int m_pow2 = m * m;
int n_pow2 = n * n;
complex<long double> down_part2 = alpha * pow(constant * (m_pow2 + n_pow2), r);
complex<long double>down_part = down_part1 + down_part2;
restored[index] = upper_part / down_part;
}
}
</syntaxhighlight>

== Последовательная сложность алгоритма ==
В пункте 1.3 мы получили, что число операций, необходимых для вычисления формулы, равно <math>15</math>. Всего формула будет применена <math>N^2</math> раз. Итого, в последовательном варианте сложность составит <math>15 * N^2</math> операций.

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

== Ресурс параллелизма алгоритма ==
Поскольку в каждой позиции <math>(i, j)</math> вычисления производятся независимо, то ресурс параллелизма будет равен размеру матриц, то есть <math>N^2</math>.

== Входные и выходные данные алгоритма ==
'''Входные данные:''' Спектр искажённого изображения и спектр ядра свёртки. Обычно дано искажённое изображение (от которого берётся преобразование Фурье). Мы считаем параметры ядра известными. Оно генерируется, и затем от него берётся преобразование Фурье (получается спектр).

'''Выходные данные:''' Спектр исходного изображения. Спектр можно преобразовать в изображение с помощью обратного преобразования Фурье.
== Свойства алгоритма ==
= ЧАСТЬ. Программная реализация алгоритма =
== Особенности реализации последовательного алгоритма ==
== Локальность данных и вычислений ==
== Возможные способы и особенности параллельной реализации алгоритма ==
== Масштабируемость алгоритма и его реализации ==

== Динамические характеристики и эффективность реализации алгоритма ==
== Выводы для классов архитектур ==
== Существующие реализации алгоритма ==

= Литература =
* Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. – М.: Наука, 1979.
* Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. – М.: Наука, 1990.
* Гудмен Дж. Введение в фурье‐оптику. – М.: Мир, 1970. 364 с.

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-05T15:09:49Z

Kiseliov: /* Схема реализации последовательного алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-05T14:10:37Z

Kiseliov: /* Схема реализации последовательного алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-05T10:11:22Z

Kiseliov: /* Схема реализации последовательного алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-05T07:51:31Z

Kiseliov: /* Математическое описание алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-04T17:07:55Z

Kiseliov: /* Математическое описание алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-04T17:06:33Z

Kiseliov: /* Ресурс параллелизма алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-04T17:06:16Z

Kiseliov: /* Последовательная сложность алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-04T17:05:50Z

Kiseliov: /* Математическое описание алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-04T17:02:46Z

Kiseliov: /* Математическое описание алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-04T17:02:02Z

Kiseliov: /* Ресурс параллелизма алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-04T16:57:26Z

Kiseliov: /* Вычислительное ядро алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-04T16:56:20Z

Kiseliov: /* Вычислительное ядро алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-04T16:56:10Z

Kiseliov: /* Последовательная сложность алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-04T16:49:16Z

Kiseliov: /* Макроструктура алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-04T16:44:18Z

Kiseliov: /* Вычислительное ядро алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-04T16:43:58Z

Kiseliov: /* Вычислительное ядро алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-04T16:42:29Z

Kiseliov: /* Математическое описание алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-04T16:42:07Z

Kiseliov: /* Вычислительное ядро алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-04T16:33:36Z

Kiseliov: /* Математическое описание алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-04T16:32:16Z

Kiseliov: /* Математическое описание алгоритма */

Участник:Kiseliov/Метод регуляризации Тихонова

2022-10-04T16:22:02Z

Kiseliov: /* Математическое описание алгоритма */