CC BY
36
2008 ВЕСТНИК САНКТ-ПЕТЕРБУРГСКОГО УНИВЕРСИТЕТА. Сер. 10. Вып. 2
ПРИКЛАДНАЯ МАТЕМАТИКА
УДК 519.6
Ю. К. Демьянович, О. М. Косогоров
О ПАРАЛЛЕЛЬНОМ ВЭЙВЛЕТНО-СПЛАЙНОВОМ СЖАТИИ НА НЕРАВНОМЕРНОЙ СЕТКЕ *)
Введение. Теория вэйвлетов (всплесков) появилась сравнительно недавно (несколько десятилетий тому назад); к настоящему времени она завоевала прочные позиции в математике (см. [1-4] и др.) и нашла глубокие приложения в физике, астрономии, медицине и, конечно, в инженерном деле, поскольку основной результат этой теории -эффективные алгоритмы обработки больших потоков информации. Под эффективностью в данном случае понимают экономное (с точки зрения экономии ресурсов компьютера: памяти и времени обработки) разложение потока информации на составляющие, среди которых выделяются основной информационный поток, уточняющий информационный поток и информационный поток с несущественной информацией. Как правило, основной информационный поток значительно менее плотный, чем исходный поток информации, потому его можно передать быстро, и при этом не требуется использовать линии связи большой ширины. Уточняющий информационный поток не во всех случаях необходим, его можно передавать фрагментарно в зависимости от потребностей. Наконец, поток с несущественной информацией вообще может быть отброшен. Конечно, вопрос о том, какая информация является основной, какая уточняющей, а какая - несущественной, выходит за рамки математических исследований и должен решаться в каждом отдельном случае специалистом предметной области.
Роль теории вэйвлетов (всплесков) состоит в том, что она дает предметному специалисту достаточно широкий арсенал средств, из которых он может выбрать то, которое ему подходит для обработки (для разложения на составляющие) интересующего его потока информации. Такими средствами в теории вэйвлетов являются наборы вложенных (основных) пространств функций и их представлений в виде прямой (а иногда и ортогональной) суммы вэйвлетных пространств. Весьма важными являются базисы основных пространств, а также базисы вэйвлетов (всплесков); построению и изучению свойств таких базисов посвящено большое количество работ.
Значительный вклад в развитие теории вэйвлетов внесли работы Мейера, Чуй, До-беши, Маллы и др. (число публикаций в этой области, по-видимому, более тысячи). Самые известные книги по этой тематике указаны в литературе (см. [1-4]); их можно рекомендовать для углубленного изучения теории, а также ее приложений.
Известно, что на двукратно измельчающейся сетке полиномиальные 1?-сплайны образуют телескопическую систему пространств, на основе которой строятся вэйвлетные
Работа частично выполнена при финансовой поддержке Российского фонда фундаментальных исследований (гранты К* 07-01-00269 и 07-01-00451).
© Ю. К. Демьянович, О. М. Косогоров, 2008
разложения ([1-5]); для равномерной сетки это устанавливается с помощью преобразования Фурье (см., например, [3]), а для неравномерной сетки - использованием специального дробно-рационального тождества ([6-8]). Применение неравномерной сетки позволяет улучшить приближение функций без усложнения вычислений. Однако для дальнейшего улучшения приближения могут понадобиться различные степени измельчения сетки в разных частях рассматриваемого промежутка: для этого двукратное измельчение недостаточно. Особую заботу представляет вэйвлетное разложение в случае неравномерной сетки, поскольку обычно применяемое на равномерной сетке преобразование Фурье в условиях неравномерной сетки использовать затруднительно. Оказалось, что применение биортогональной системы функционалов позволяет построить вэйвлетные разложения при произвольном измельчении сетки (это ведет к упрощениям и в случае равномерной сетки).
В данной работе даны параллельные алгоритмы для вэйвлетного разложения телескопической системы пространств сплайнов второй степени при произвольном способе измельчения сетки (вывод формул см. в работе [8]). Общий ход рассуждений состоит в следующем. Первоначально с каждой сеткой связываются £?-сплайны второй степени. Далее устанавливается, что если новая сетка получена из исходной удалением одного узла, то пространство минимальных сплайнов, построенных для крупной сетки, содержится в пространстве аналогичных сплайнов, построенных для мелкой сетки. Использование биортогональной (к В-сплайнам) системы функционалов позволяет получить прямое сплайн-вэйвлетное разложение и соответствующие формулы реконструкции и декомпозиции; отсюда следует, что последовательное удаление или добавление узлов в исходную сетку приводит к телескопическим системам пространств рассматриваемых сплайнов и к прямому разложению исходного пространства. От исходной сетки можно перейти к любой более мелкой сетке добавлением любой конечной совокупности узлов, не совпадающих с узлами расширяемой сетки (между любыми двумя узлами исходной сетки число добавляемых узлов может быть различным в зависимости от выбранной пары соседних узлов). Формулы декомпозиции и реконструкции получаются при помощи соответствующих биортогональных систем функционалов. Все результаты верны в случаях конечного отрезка и бесконечного интервала.
Предлагаемые в работе алгоритмы позволяют полностью загрузить двухъядерную или двухпроцессорную архитектуру. Их целесообразно использовать для обработки больших потоков числовой информации с резко меняющимися характеристиками. В первом приближении упомянутые числовые потоки характеризуются наличием регулярной и сингулярной частей. Сжатие информации может быть осуществлено путем отбрасывания потока с несущественной информацией (при этом важно грамотно выделить поток с несущественной информацией из исходного информационного потока в каждом отдельном случае). Алгоритм восстановления позволяет провести восстановление с потерей (таким образом, сжав информацию) или (при необходимости) полностью восстановить исходный числовой поток.
Декомпозиция и реконструкция. Построение вэйвлетного разложения определяется формулами декомпозиции. Приведем их для исходного потока С = {с*}, каждое значение с* которого поставлено в соответствие узлу XI сетки X = {х{\. Результатом разложения являются поток А = {а^} и вэйвлетные коэффициенты В = {6^}. Для краткости здесь рассматривается лишь алгоритм отбрасывания одного узла Хк и соответствующей ему величины Ск', естественно, что при этом поток А существенно отличается от потока С лишь одним значением, а вэйвлетные коэффициенты все равны нулю, за исключением коэффициента Ък- Отбрасывание большего количества узлов
можно получить последовательным применением приводимого алгоритма.
4
а* = с* при / $ к — 2.
_ 2^+2 - Ж/1—1 Хк+2 - Жк+1
®>к — 1 Ск — Х Ск — 2]
Хк+1 - Хк-1 жк+1 - Ж*_1
а» = с»+1 при г ^ к,
Ь, = 0 при ; ф к,
, Ж*+3 - Жй+1 / Хк+2 - Ж*+1 Хк+2 - Ж*_1 \ , Ж*+1 - Ж*;
о* = ------------- ------------ С*_2-------------------------С*_1 + Ск-------Ск+1-
Жй+1 - Ж/1—1 \ Ж*+3 - Хк Хк+3 - Хк / Ж*+3 - Жк
Для реконструкции необходимо использовать формулы
= о,- при ] ^к — 2,
Жй+1 - Жк_1 Хк+2 - Жй+1
сд-1 —--------------ад-1 Н--------------а*—г 9
жд+г - ж*_1 ж*+2 - жк_1
Жй+1 - Хк 1 Хк+3 - Хк+1 , ,
Ск =---------- ак Н------------------а*_1 + Од,
жд+з - жк хк+з - жк
с.? = ПРИ ^ + 1-
Замечание. Формулы декомпозиции и реконструкции даны в уравновешенной форме: хотя количество арифметических операций можно уменьшить за счет вынесения некоторых сомножителей за скобки, однако такое преобразование может привести к неустойчивым вычислениям; в данном случае количество операций возрастает незначительно, но резко повышается устойчивость счета.
Распараллеливание вычислений. Рассмотрим распараллеливание вычислений при построении вэйвлетного разложения и при восстановлении исходного потока. Представим порядок вычислений при использовании двух параллельных процессов. Распараллеливание демонстрируется с помощью параллельных форм, в которых ярусы нумеруются последовательно; на каждом ярусе указываются: две очередные вычислительные операции, используемая память (в числе единиц хранимых данных) и число необходимых транзакций (указывается сумма числа пересылок между процессами и числа пересылок при поступлении данных).
Замечание. В случае применения равномерной сетки формулы значительно упрощаются, однако такую сетку целесообразно использовать лишь на начальном этапе сжатия, когда числовой поток может быть идентифицирован с равномерной сеткой. В процессе обработки появляется сильно разреженная неравномерная сетка, и потому дальнейшее сжатие проводится по приведенному выше параллельному алгоритму. Программная реализация декомпозиции и реконструкции показала быструю и эффективную работу алгоритма на ряде модельных примеров.
Теперь поясним данные, приведенные в табл. 1. Пусть для каждого процессора выделено определенное количество ячеек памяти. Будем считать, что на каждом ярусе
Ярус I процесс II процесс Пам Тр
1 Хк+2 - £к+1 ^Л+З - хк 6 4
2 Хк+2 - %к-1 Хк+1 - Хк-1 8 3
3 хк+2-хк+1 хк+3-хк Хк+2 - Хк-1 Хк+з - хк 8 2
4 Хк+2 - £Л + 1 «Л = Ск- 2 Хк+з - хк %к+2 - Хк-1 '»к = Ск-1 Хк+3 - %к 11 2
5 «Л - чк %к+3 - %к+1 9 1
6 Хк+1 - хк %к+3 - %к+1 Хк+1 - Хк-1 7 1
7 Хк+1 - хк Хк+з - хк %к+3 - %к+1 ■Шк = ик V,. Хк+1 - хк-1 5 1
8 хк+1 - (-'к-1“ 1 ^Л+З - хк ск +ык 7 2
9 %к+2 - Хк-1 Хк+1 - Хк-1 Хк+2 - £к + 1 + 1 - Хк-1 4 2
10 %к + 2 - Хк-1 Рк = Ск-1 Хк+1 - хк-1 Хк+2 - Хк+1 <1к — Ск-2 Хк+1 - хк-1 3 2
11 о-к-1 =Рк-Як , , | . хк+1 - хк Ок = (Ск+и>к) Ск + 1 Хк+з - хк 2 2
(итерации) соблюдается следующий порядок действий: 1) поступление данных на процессор; 2) запись в память; 3) пересылки с одного процессора на другой; 4) стирание из памяти; 5) запись результатов вычислений с процессора в память. Также положим, что если результат вычислений на процессоре нужен только для вычислений на этом же процессоре на следующей итерации, то он может не записываться в память, т. е. будем считать данную информацию «уже поступившей» на процессор.
Для наглядности представим состояние двухпроцессорной системы на каждом п-м ярусе параллельных вычислений в виде табл. 2.
Таблица 2. Состояние системы на ярусе п
Состояние I процесс II процесс
Поступление данных сразу на процессор 1
Поступление данных с записью в память 2
Пересылки между процессорами 3
Стираемые данные 4
Данные с процессора 5
Свободные ячейки 6
В клетке табл. 2, отмеченной цифрой 1, указываются данные, поступающие извне для вычислений на первом процессоре. Они не используются далее, поэтому не записываются в память, а поступают непосредственно на процессор. В клетке 2 приводятся данные, поступающие извне для вычислений на первом процессоре, которые понадобятся в дальнейшем, поэтому они записываются в память, отведенную для первого процессора (в скобках указан номер ячейки памяти, в которую производится запись). В клетке 3 даются данные, пересылаемые со второго процессора (или из памяти, отведенной для него) на первый (если эти данные записываются в память, то в скобках указывается номер ячейки памяти; если же они не понадобятся для вычислений на
следующих итерациях, то в скобках написано «на пр»), В клетке 4 приводятся данные, которые более не нужны для дальнейших вычислений и стираются из памяти (в скобках дан номер освобождаемой ячейки памяти). В клетке 5 указывается результат вычислений на первом процессоре на данной итерации. Он может записываться в память (тогда в скобках приводится номер занимаемой ячейки памяти); если же результат вычислений на первом процессоре нужен только для вычислений на нем на следующей итерации, то он в память не записывается (тогда в скобках написано «на пр»), В клетке 6 указываются номера ячеек памяти, свободных перед следующей итерацией. Все вышесказанное аналогично для второго процесса.
В рассматриваемом случае на каждом ярусе нам понадобится не более шести ячеек памяти, отведенных для каждого из двух процессоров. Получим такое состояние системы:
на ярусе 1:
на пр
с зап Хк+2 (1) £Л+1 (2) жЛ+3 (1) хк (2)
перес
стир
с проц Хк+2 - £Л+1 (3) ^к+з - хк (3)
своб яч 4 5 6 4 5 6
на ярусе 2:
на пр Хк- 1 Хк- 1
с зап
перес £к+1 (4)
стир Хк+2 (1)
с проц Хк+2 - %к-1 (1) £к+1 - хк-1 (5)
своб яч 4 5 6 6
на ярусе 3:
на пр
с зап
перес ^к+з - хк (4) Хк+2 - хк-1 (на пр)
стир ^к+з - хк (3)
с проц Хк+2 - £к+1 , . (на пр) ^к+з - хк Хк+2 - Хк-1 , . (на пр) Хк+3 - Хк
своб яч 5 6 3 6
на ярусе 4:
на пр
с зап с-к-2 (5) Ск-1 (3)
перес
стир
с проц Хк+2 - £Л+1 , . ик = ск-2 (на пр) ^Л+З - Хк Хк-\-2 Хк—1 //->\ Ък = Ск-1 (6) Хк+з - %к
своб яч 6 нет
на пр
с зап
перес Хк±2-%к-1 , v vk = cfc_ 1 (на пр) хк+3 - хк
стир жк+3 (1); жк+1 (4); vk (6)
с проц «к - Vk (6) ^fc+з - xk+i (на пр)
своб яч Нет 1 4 6
на ярусе 6:
на пр
с зап
перес хк (на пр)
стир £fc+i (2) хк (2)
с проц £fc+i - хк (на пр) £fe+3 - Хк+1 , ч (на пр) хк+1 - Хк_1
своб яч 2 12 4 6
на ярусе 7:
на пр
с зап
перес ик - vk (на пр)
стир Хк+з - Хк (4); ик - vk (6)
с проц Хк+1 ~Хк , ч (на пр) Zfc+З - хк ™к = ик vk (на пр) хк+і - хк-г
своб яч 2 4 6 12 4 6
на ярусе 8:
на пр Cfe + 1 Cfc
с зап
перес
стир
с проц Хк+1 - хк „ Cfc+1 (2) Хк+3 - хк cfc + wk (1)
своб яч 4 6 2 4 6
на ярусе 9:
на пр
с зап
перес £fc+i - хк-1 (на пр) Хк+2 - £fc+i (на пр)
стир Хк+2 - Хк-1 (1); Хк+2 - Хк + 1 (3) Xk+l - xk-i (5)
с проц Хк+2 - £fc_l , . (на пр) Xk+i - xk-i Хк+2 - £fc+i , . (на пр) хк+1 - Хк-1
своб яч 13 4 6 2 4 5 6
на ярусе 10:
на пр
с зап
перес ск-1 (на пр) Ск-2 (на пр)
стир Cfc-2 (5) Ск-1 (3)
с проц Xk+2-Xk-l Рк = ск-1 (на пр) Xk+i - Хк-г Хк+2 ^fe+І Як = ск-2 (2) Xk+i - Хк-г
своб яч 1 3 4 5 6 3 4 5 6
на пр
с зап
перес £fe+2-a;fe+1 Qk = Cfe-2 (на пр) хк+г - хк_1 Xk+1-Xk , . ск+1 (на пр) £fc+3 - хк
стир хк±1-хк ск+1 {*) хк+3 - хк Ск +wk (1) qк (2)
с проц о-к-1 =Рк-Як (1) г. / 1 \ Хк ,л \ Ok = (Cfc + wk) Ск + 1 (1) Zfc+з - Хк
своб яч 2 3 4 5 6 2 3 4 5 6
Рассмотрим также возможный вариант распараллеливания вычислений при реконструкции исходного потока. Данные табл. 3 можно пояснить аналогичным способом (с помощью табл. 2).
Таблица 3. Распараллеливание формул реконструкции
Ярус I процесс II процесс Пам Тр
1 £fc+i - Хк- 1 Хк+з - Хк + 1 6 4
2 Хк+2 - Хк-1 Хк+з - Хк 10 2
3 £fc+i - Хк-1 Хк+2 - Хк-1 Хк+з - Хк + 1 Хк+з - Хк 6 0
4 £fc+i - хк-1 «ft = ttfc-1 Хк+2 - £fc_l Хк+з - Хк + 1 Wk = О.к-1 Хк+з - Хк 7 2
5 Хк+2 - £fc + l Pk = wk + Ьк 6 1
6 Хк + 2 - £fc + 1 Хк+2 - Хк-1 Хк+1 - Хк 3 0
7 Хк + 2 - £fc + l «ft = a-к-2 Хк+2 - Хк-1 Хк+1 - Хк Хк+з - Хк 2 1
8 Ск-1 =tlk+t>k Хк+1 - Хк Як = ак Хк+з - Хк 2 1
9 Ск=Рк + Як 2 0
Замечание. В приведенных выше параллельных формах ставилась задача минимизировать количество транзакций (в первую очередь при пересылке данных между процессорами) при минимально возможном количестве ярусов вычислений.
Summary
Demjanovich Yu. К., Kosogorov О. М. On parallel spline-wavelet compression on irregular grid.
Parallel forms of decomposition and reconstruction algorithms in case of quadric -B-splines on irregular grid are presented. These presented parallel forms make it possible to use a two-processor system fully.
Литература
1. Петухов А. П. Введение в теорию базисов всплесков. СПб.: Изд-во СПбГТУ, 1999. 132 с.
2. Добеши И. Десять лекций по вейвлетам / Пер. с англ. Е. В. Мищенко; Под ред. А. П. Петухова. Москва; Ижевск: НИЦ «Регулярная и хаотическая динамика», 2004. 464 с.
3. Чуй К. Введение в вэйвлеты / Пер. с англ. Я. М. Жилейкина. М.: Мир, 2001. 412 с.
4. Малла С. Вэйвлеты в обработке сигналов / Пер. с англ. Я. М. Жилейкина. М.: Мир, 2005. 671 с.
5. Новиков И. Я., Стен кин С. Б. Основы теории всплесков // Успехи матем. наук. 1998. Т. 53, № 6 (324). С. 53-128.
6. Демьянович Ю. К. Калибровочное соотношение для В-сплайнов на неравномерной сетке // Матем. моделирование. 2001. Т. 13, № 9. С. 98-100.
7. Демьянович Ю. К. Всплесковые разложения в пространствах сплайнов на неравномер-
ной сетке // Докл. РАН. 2002. Т. 382, № 3. С. 313-316.
8. Демьянович Ю. К. Локальный базис всплесков на неравномерной сетке // Зап. науч.
семинаров Петерб. отд. Матем. ин-та РАН. 2006. Т. 334. С. 84-110.
Статья рекомендована к печати проф. А. Н. Тереховым.
Статья принята к печати 4 декабря 2007 г.
