CC BY
53
ВЕСТНИК ТОМСКОГО ГОСУДАРСТВЕННОГО УНИВЕРСИТЕТА 2014 Управление, вычислительная техника и информатика № 4 (29)
УДК 519.6
Б.М. Шумилов, У.С. Ыманов
СПЛАЙН-ВЕЙВЛЕТЫ, ОРТОГОНАЛЬНЫЕ МНОГОЧЛЕНАМ, И ОПТИМИЗАЦИЯ ВЫЧИСЛЕНИЙ ВЕЙВЛЕТ-ПРЕОБРАЗОВАНИЯ
Исследование выполнено при финансовой поддержке РФФИ в рамках научного проекта № 14-31-50353-мол_нр.
Для случая эрмитовых сплайнов произвольной нечетной степени исследовано построение системы базисных
сплайн-мультивейвлетов, реализующих условия ортогональности многочленам той же степени. Рассмотрены
варианты построения и обращения блока фильтров. Приведены результаты численных экспериментов.
Ключевые слова: мультивейвлеты нечетной степени; ортогональность многочленам.
Кубические сплайны гладкости С2 завоевали общее признание в среде инженеров-дорожников как адекватный способ математического представления пикетного метода трассирования реконструируемых автомобильных дорог [1-4]. Перспективы использования здесь эрмитовых сплайнов состоят в том, что они позволяют в явном виде, через значения сплайн-коэффициентов, учесть геометрические ограничения как на контрольные точки (пикеты автомобильной трассы), так и на тангенсы трассы (заданные направления касательных при въезде на мостовое сооружение или примыкание) и радиусы закруглений (заданные значения кривизны на разгонных и тормозных участках трассы). В работах [5, 6] были построены мультивейвлеты эрмитовых сплайнов пятой и третьей степени, которые ортогональны любым многочленам той же степени. В представленной статье получено матричное решение для коэффициентов мультивейвлетов произвольной нечетной степени, ортогональных многочленам той же степени. Предложена эффективная его модификация, которая позволяет воспользоваться для численного решения хорошо разработанным методом матричной прогонки [7. С. 103].
1. Построение системы базисных сплайн-мультивейвлетов на конечном отрезке
Пусть на отрезке [а, Ь] задана вложенная последовательность равномерных сеток АЬ : хI = а + И • /, / = 0, 1,..., 2Ь, И = (Ь - а) / 2Ь, Ь > 0. Если базисные функции КЬ/кк (х) = фк (V - /), к = 0, 1, ..., г V/, где V = (х - а) / И с центрами в узлах сетки АЬ порождены сжатиями и сдвигами г+1 функций вида [8. С. 82]:
Фк (Х) = ^("1)к Юк (-')при-1 < х < а
(х) при 0 < X < 1,
г+1 г-к (г + В)1
где &к (х) = (1 -X)г+ -— хк+р, к = 0,1,...,г, то, при условии отсечения выступающих за концы от-
р=0 к!Р!г!
резка половинок функций ф0(х), ..., фг(х), полученное пространство УЬ является пространством эрмитовых сплайнов степени 2г+1 гладкости С.
Поскольку сетка АЬ-1, Ь > 1, получена из АЬ посредством удаления каждого второго узла, то соответствующее пространство УЬ-1 с базисными функциями КЬ-11,к с носителями, в два раза большими по ширине, и центрами в четных узлах сетки АЬ вложено в УЬ. Остаток от разности пространств УЬ и УЬ—1 размерности (г+1) (2Ь+1-2Ь-1-1) = (г+1)2Ь-1 называется пространством мультивейвлетов. Будем искать базисные мультивейвлеты МЬ/к (х), к = 0, 1, ..., г V/, как линейные комбинации базисных эрмитовых сплайнов на сетке АЬ+:, удовлетворяющие условиям ортогональности многочленам 2г+2-го порядка, т. е.
b
J Mfk (x)xmdx = 0, k = 0,1,..., rVi (m = 0,1,...,2r +1), (1)
a
и имеющие минимальные из возможных носители.
Теорема. Пусть задан отрезок [0, 2L+1], L > 0, с сеткой AL+1 : x, = i, i = 0, 1,..., 2L+1, с единичным шагом h = 1, и разложения базисных мультивейвлетов имеют вид
MLk (x) = £ £ ajФг (2t - j), t = (x - x2i), -1 < t < p +1. (2)
l=0 j=0
Тогда при условии p = 2 три матрицы R0, R1, R2 размерности (2r + 2) x (r +1) заданы элементами
j+1
Rj,i = J Ф1 (2t- j)tmdt, j = 0,1,2, l = 0,1,.,r, m = 0,1,...,2r +1, j-1
каждый из r+1 столбцов блочной матрицы |A0nner /A'2iner J = -[R0 | R2]-1 R1 дает значения коэффициентов aj, j = 0,2, l = 0,1,.,r, соответствующего k-го базисного мультивейвлета, полностью лежащего внутри отрезка [0, 2L+1], L > 1, а коэффициенты aj = {1, l = k; 0, l Ф k}. При L = 0 элементы матриц R0, R2 вычисляются в пределах интервала [0, 2], а матрица коэффициентов aj, j = 0,2, l = 0,1,., r, обозначены
„ Г л center / a center "I г F)0 i n2i-1 г)1
как | A 0 /A 2 1 = -[ R | R ] R .
При L > 0 для крайних слева базисных мультивейвлетов элементы матрицы R0 вычисляются по укороченному слева интервалу
1
Rm,l = J Ф1 (21)tmdt, l = 0,1,., r, m = 0,1,. ,2r +1, 0
а коэффициенты разложения (2) aj, j = 1,2, l = 0,1,.,r, даются значениями столбцов матрицы ^ A}eft /A'2eft J = -[R1 | R2 ]-1 R0 при условии, что коэффициенты a0 = {1, l = k; 0, l Ф k}. Для крайних справа базисных мультивейвлетов элементы матрицы R2 вычисляются по укороченному справа интервалу
2
Rm,l = J Ф1 (2t - 2)tmdt, l = 0,1,., r, m = 0,1,.,2r +1, 1
а коэффициенты разложения (2) aj, j = 0,1, l = 0,1,.,r, даются значениями столбцов матрицы [A0ighVA[ight J = -[R0 | R1]-1 R2 при условии, что коэффициенты a2 = {1, l = k; 0, l Ф k}.
Система функций MLik (x), k = 0, 1,...,r, i = 1,2,.,2L, удовлетворяет условиям (1) с носителями не более чем из двух шагов сетки AL+1 и образует базис в пространстве WL, L > 0.
Доказательство. а) Пусть L > 1 и носители мультивейвлетов [x2i-i, x2i+p+i], p > 1, полностью располагаются внутри отрезка [0, 2L1]. Приp > 1 для всех i = 1,2,.,2L-1-[p/2] имеет место разложение (2). Подставим его в (1) и вычислим все необходимые интегралы, учитывая при этом, что подынтегральные выражения обращаются в нуль вне носителей вейвлетов. Условия ортогональности (1) определяют однородную систему 2r+2 уравнений метода моментов [9. С. 42] относительно коэффициентов мультивейвлета MLik(x). В силу линейной независимости на интервалах [x2i-1, x2i+p+1] пробных функций
N^ k (x), l = 0,., r, j = 0,1,., p, и поверочных функций xm, m = 0,1,., 2r+1, ранг полученной системы равняется min[2r+2, (r+1)(p+1)]. Если предположить, что носитель мультивейвлета равен трем шагам
сетки Д1+1, т.е. р = 1, то однородная система становится квадратной и, будучи невырожденной, может иметь только тривиальное решение. Поэтому будем далее считать, что носитель мультивейвлета равен четырем шагам сетки Д1+1, т.е. мультивейвлет построен из 3г+3 базисных сплайнов. В этом случае ранг матрицы равен 2г+2. Поэтому однородная система 2г+2 уравнений с 3г+3 неизвестными имеет г+1 линейно независимых частных решений. При этом количество полученных для каждого номера 1 мульти-вейвлетов, полностью лежащих внутри отрезка [0, 2м], равно (г+1)(2ь-2) = (г+1)2ь- 2(г+1), что на 2(г+1) меньше разности между размерностями пространств Уь+1 и Уь\ (г+1)(21+1+1) --(г+1)(2^+1) = (г+1)2ь.
Вблизи концов отрезка интервалы интегрирования не выходят за пределы отрезка. Поэтому с учетом того, что по краям отрезка от выступающих за концы отрезка функций ф0 (¿), ..., фг (¿) остается по половинке, граничные мультивейвлеты отличаются от мультивейвлетов, предложенных выше, при условии ортогональности многочленам 2г+1-й степени. Аналогичные рассуждения приводят к двум системам линейно независимых функций М1^ (х), к = 0, 1,...,г, отдельно на левом и на правом концах, дополняющих до базиса построенную выше мультивейвлет-систему.
б) Пусть 1 = 0. Вычисления снова дают г+1 линейно независимых частных решений, которые образуют базис в пространстве У - У0 с размерностью (г+1)3-(г+1)2 = (г+1).
в) Окончательно, в соответствии с правилами линейной алгебры выписываются матричные формулы для вычисления коэффициентов эрмитовых мультивейвлетов нечетной степени. Теорема доказана.
2. Построение и обращение блока фильтров
Если упорядочить базисные сплайн-функции в виде единой матрицы-строки
ФЬ =
<0, <Д,..., Кг, <0, <1,..., N21
2Ь ,г
, то можно записать функции фЬ 1 в виде линейных комбинаций функций ф : ф - = ф
: р.
Здесь блоки матрицы составлены из коэффициентов масштабных соотношений (также для единичного шага) [10]:
ф0(^) ф0(2'- -к)
фЛ() = Х нк к=0 фг(2Г- -к)
_фг С) _ _фг(2Г- -к)_
где Н2 = и 1Ли, И = diag(1,2 1,...,2 г), Н0 = БН2Б 1 и матрица и размерности (г +1) х (г +1) задана
элементами
, чг+1+к _ ,■ (г +1 + к)! и к, ] =(-1)г+1+к V.,., , к, ] = 0,1,., г,
(г +1 + к - ]) а Л = ^ (2-г-1,...,2-2г-1), Б = ^ (1,-1,...,(-1)-г).
Аналогично запишем базисные мультивейвлет-функции 2г+1-й степени на уровне разрешения
Ь в виде матрицы-строки у1 =
М^, МЬЛ,..., ,..., М Ьг
. Тогда для уровня разрешения Ь-1
можно выразить функции у1 1 в виде линейных комбинаций функций фЬ: у1 1 = фЬ $, где блоки матрицы $ составлены из коэффициентов разложений (2). Соответствующие коэффициенты
т
сплайна будем собирать в вектор с =
1,0 ^.1,1
рЬ,0 р
1,г
-0 , С0 , С0 , С1 ,..., С2Ь
, и соответствующие муль-
тивейвлет-коэффициенты - в вектор Б1 = Д1,0, Д1,1
,Д1,г,..., ДТ
. Тогда с использованием
обозначений для блочных матриц процесс получения С1 из С1 1 и Д 1 может быть записан так [11. С. 101]:
СЬ = [рЬ
]
С
Ь-1
П
Ь-1
(3)
Разрешимость системы (3) относительно С - , гарантирована линейной независимостью базисных функций. Для упрощения вычислений матрицу [Р | $] предлагается сделать блочной трехдиа-гональной, изменив порядок неизвестных так, чтобы блоки матриц Р и $ перемежались:
к V-1 = сЬ,
иЬ =
С0 ,
С,
Ь,1
> С0 , П1 , П1 ,
п1 , С1 ,
пь ,г с > 21 '
Ь,г 2ь
к1 =
Я!
Я 2Т о
л сеМег Л0
о
Я 0Т
А2еп1ег Я1
кЬ =
Г Я I 0
Я 2Т 41ей Я 0Т
0 А* Я\
0 0 Я Т
0 0 0
0 0
л тпег А0
А2
Я
0
0
0
0
0 0
0 I
0 ^ 0 0
Я 0 Я
Ь > 1.
1 У
Здесь О обозначает матрицу г+1-го порядка с нулевыми коэффициентами, тогда как I - единичная матрица г+1-го порядка, многоточия обозначают повторяющиеся блоки матрицы КЬ. При этом с целью компенсации единичного шага сетки в уравнениях (3) в качестве исходных СЬ нужно использовать значения функции и производных, домноженные на Н в соответствующей степени: {/к)(гН)Нк, к = 0, 1, ..., г, / = 0, 1,..., 2Ь}, всего (г+1)(2Ь+1) чисел.
Процедуру разбиения СЬ на часть соответствующую низшему разрешению, и уточняющие коэффициенты ПЬ-1 можно применить рекурсивно и к самой части СЬ-1. Следовательно, исходные значения СЬ можно представить в виде иерархии все более грубых версий с разрешениями С0, С1, ..., СЬ-1 и уточняющих деталей П0, П1, ..., ПЬ-1. При этом по величине вейвлет-коэффициентов П], у = 0, 1,., Ь-1, можно судить об их значимости для результирующей сплайн-аппроксимации. Незначимые убираются с целью сжатия информации.
3. Пример сравнения мультивейвлетов степени 7 и 11 при одном и том же числе базисных функций 12
Для хе [0, 1] рассмотрим в качестве примера приближение функции Хартена [12]:
^т(3л:х), х < -3,
" 1 < х <
I ( х ) =
|8т(4лх)|, з < х < -3,
- -2 8т(3л:х),
х > -3.
Если для случая г = 5 аннулировать все вейвлет-коэффициенты, оставив только сплайн-коэффициенты С0 = [0,2968; -32,14; 2149; -7,73 1 04; 1,854 106; -3,122 107; 0,3851; 55,15; 3146; 1,222 1 05; 3,31 106; 6,49 1 07] т, то получится некоторый сглаживающий многочлен 11-й степени, весьма незначительно отличающийся от МНК-решения (рис. 1), представленного коэффициентами [0,1818; -22,93; 1694; -6,193 • 104; 1,485-106; -2,48 107; 0,2333; 42,72; 2514; 9,995 1 04; 2,74-106; 5,407- 107]Т.
Рис. 1. Сравнение значений сглаживающего многочлена 11-й степени (точки) и МНК-многочлена 11-й степени (сплошная линия), построенных по значениям функции Хартена (жирная линия)
Для случая г = 3 остаются только сплайн- коэффициенты С1 = [0,359; -16,22; 448,8; -6532; 0,293; -1,104; 53,74; 570,8; 0,57; 32,56; 755,3; 1,137 104]т. В результате получается сглаживающий сплайн 7-й степени с одним узлом посередине отрезка. Отличие от МНК-решения с коэффициентами [0,2312; -12,22; 377,1; -5711; 0,3952; -0,2757; 38,69; 229,2; 0,2653; 19,64; 463,6; 7208]т выглядит более значительным при заметном улучшении подгонки в окрестности излома функции (рис. 2).
Рис. 2. Сравнение значений сглаживающего сплайна 7-й степени (точки) и МНК-сплайна 7-й степени (сплошная линия), построенных по значениям функции Хартена (жирная линия)
Заключение
В работе представлена общая схема построения эрмитовых сплайн-вейвлетов, ортогональных многочленам. Полученные результаты предоставляют широкие возможности для оптимизации методов обработки численной информации, которые допускают и параллельную реализацию [13. С. 139; 14].
ЛИТЕРАТУРА
1. Бойков В.Н., Шумилов Б.М. Сплайны в трассировании автомобильных дорог. Томск : Изд-во ГУ Томский ЦНТИ, 2001. 164 с.
2. Система проектирования МогСАБ. Проектирование автомобильных дорог: руководство пользователя / И.В. Кривых,
Д.А. Петренко, В.Н. Бойков и др. 2-е изд., испр. Томск : Изд-во Том. ун-та, 2010. 250 с.
3. БойковВ.Н. САПР автодорог - перспективы развития // САПР и ГИС автомобильных дорог. 2013. № 1(1). С. 6-9.
4. Петренко Д.А. Новое поколение программных продуктов ИндорСофт // САПР и ГИС автомобильных дорог. 2013. № 1(1).
С. 10-17.
5. Шумилов Б.М., Кудуев А.Ж. Новый тип мультивейвлетов пятой степени, ортогональных многочленам пятой степени // Вест-
ник Томского государственного университета. Управление, вычислительная техника и информатика. 2012. № 4(21). C. 108116.
6. Шумилов Б.М. Мультивейвлеты эрмитовых сплайнов третьей степени, ортогональные кубическим многочленам // Матема-
тическое моделирование. 2013. Т. 25, № 4. С. 17-28.
7. Самарский А.А., НиколаевЕ.С. Методы решения сеточных уравнений : учеб. пособие. М. : Наука, 1978. 591 с.
8. Завьялов Ю.С., Квасов Б.И., Мирошниченко В.Л. Методы сплайн-функций. М. : Наука, 1980. 352 с.
9. Флетчер К. Численные методы на основе метода Галёркина : пер. с англ. М. : Мир, 1988. 352 с.
10. Strang G., Strela V. Short wavelets and matrix dilation equations // IEEE Trans. Signal Processing. 1995. V. 43, Nfo. 1. P. 108-115.
11. Столниц Э., ДеРоуз Т., Салезин Д. Вейвлеты в компьютерной графике : пер. с англ. Ижевск : Регулярная и хаотическая динамика, 2002. 272 с.
12. Arandiga F., Baeza A., Donat R. Discrete multiresolution based on hermite interpolation: computing derivatives // Communications in Nonlinear Science and Numerical Simulation. 2004. №э. 9. P. 263-273.
13. Ортега Дж. Введение в параллельные и векторные методы решения линейных систем. М. : Мир, 1991. 367 с.
14. Высокопроизводительные вычисления на кластерах : учеб. пособие / под ред. А.В. Старченко. Томск : Изд-во Том. ун-та, 2008. 198 с.
Шумилов Борис Михайлович, д-р физ.-мат. наук, профессор. Е-таД: [email protected] Томский государственный архитектурно-строительный университет Ыманов Улукбек Сайдакрамович. Е-таП: [email protected]
Ошский государственный университет Поступила в редакцию 8 сентября 2014 г.
ShumilovBorisM., Ymanov UlukbekS. (Tomsk State University. Russian Federation, Osh State University. Kyrgyzstan). Spline-wavelets, orthogonal to polynomials, and optimization of calculations of wavelet-transformation. Keywords: Multiwavelets of odd degree, orthogonality to polynomials.
r 2L
For the space of the Hermitian splines of any odd degree 2r+1 of a kind SL (x) = ^ ^ CL,k NLk (x), a < x < b, with a uni-
k=0 i=0
form grid of nodes AL : Ui = a + (b - a) i / 2L, i = 0, 1,..., 2L, L > 0, and the basic functions NLk(l) (u} ) = 5j • 5k, l = 0,1,..., r, with the centers in integers, it is proposed to use as wavelets the functions MLik(x), satisfying the conditions of orthogonality to all polynomials of the 2r+2-nd order, i.e. f ML'k(x)xmdx = 0, k = 0,1,.,r Vi (m = 0,1,...,2r +1). For wavelets with the centers in even integers
J a
and the supports, which are equal to the supports of basic splines on a grid AL, the formulas for calculating the coefficients
..., D
,L-1,r
..., D
L-1, r
in the thinned grid ÀL 1 from the spline coefficients С0¡'O, С0¡'í,..., С0¡'r, ^
-rL,1
-rL,r
-rL,0
..., С , in a
' 2l
dense grid AL in the form of the solution of the linear algebraic equations system with a block three-diagonal matrix of a kind
k v-1 = cl ,
,,L _ |>L,0 сL,1 U _|_С0 , С0 ,
Г Hl
KL _
/~*L,r r\L,0 pjL,1 ' С0 , D1 , D1 ,
, D1L,r, C1L,0,
A
" H1 д center Л0 O '
k1 _ HI I h O2
O center Л 2 H1
I
left
0
h2 Aj H
left 2
0 0
Hl
H T2
AO
0
I
DL
0
right
0 ^
0
0
0
h2 Aright HOT
H1
■ CLLr ]J
L >1.
are received. Here H2 = U 1AU, H1 = diag(l,2 1,...,2 r), H0 = SH2S and the matrix Uof dimensions (r +1) x (r +1) is set
by the elements U, , =(- 1\r+i+k~J f (Г +1 + k)! , k, j = 0,1,...,r, and Л = diag(2"r-1,...,2-2r-1),
k,j \ f (r +1 + k - j)! 6V '
S = diag (1, —1,..., (—1) r ) , O designates a matrix of the r+1-st order with zero coefficients whereas I is a single matrix of the r+1-st
order, dots designate the repeating biocks of the matrix kl. The biocks of a matrix [aO™
er /a2nner] = -[R0|R2]-1 R1 are set by the
elements, respectively
j+1
RJmJ = J 4>l(2t-j)tmdt, j = 0,1,2,1 = 0,1,...,r, m = 0,1,.„,2r +1,
j-1
where
фк (t ) = 1)k Юк (- )пРи-1 ^ t ^ 0,
[ rak (t) при 0 < t < 1,
and юк (t) = (1 -1)r+1 У tk+P, k = 0,1,., r.
kW V ' ptO k ipi r!
Results of numerical experiments for different degree of a spline at the same number of basic functions in comparison with a classical method of the least squares (MLS) are presented. The improvement of adjustment on condition of coincidence of a spline knot to a break of the approximated function is revealed. It is shown that for lack of spline knots the MLS-solution and the wavelet-decision are close to each other. In the presence of spline knot the MLS -solution and the wavelet-decision differ from each other more considerably.
REFERENCES
1. Boykov V.N., Shumilov B.M. Splayny v trassirovanii avtomobil'nykh dorog [Splines in tracing of highways]. Tomsk: State Admin-
istration of Tomsk Centre of Scientific Research Publ., 2001. 164 p.
2. Krivykh I.V., Petrenko D.A., Boykov V.N. et al. Sistema proektirovaniya IndorCAD. Proektirovanie avtomobil'nykh dorog: rukovod-
stvo pol'zovatelya [System of design IndorCAD. Design of highways: user's guide]. Tomsk: Tomsk State University Publ., 2010. 250 p.
3. Boykov V.N. CAD systems for roads - evolution prospects. SAPR i GIS avtomobil'nykh dorog - CAD and GIS for roads, 2013,
no. 1(1), pp. 6-9. (In Russian).
4. Petrenko D.A. The new generation of IndorSoft software products. SAPR i GIS avtomobil'nykh dorog - CAD and GIS for roads,
2013, no. 1(1), pp. 10-17. (In Russian).
5. Shumilov B.M., Kuduev A.Zh. New type multiwavelets of the fifth degree orthogonal to quintic polynomials. Vestnik Tomskogo
gosudarstvennogo universiteta. Upravlenie, vychislitel'naya tekhnika i informatika - Tomsk State University Journal of Control and Computer Science, 2012, no. 4(21), pp. 108-116. (In Russian).
6. Shumilov B.M. Multiwavelets of the third-degree Hermitian splines orthogonal to cubic polynomials. Matematicheskoe modeliro-
vanie - Mathematical Models and Computer Simulations, 2013, vol. 25, no. 4, pp. 17-28. DOI: 10.1134/S2070048213060100
7. Samarskiy A.A., Nikolaev E.S. Metody resheniya setochnykh uravneniy [Methods of solving grid equations]. Moscow: Nauka Publ.,
1978. 591 p.
8. Zav'yalov Yu.S., Kvasov B.I., Miroshnichenko V.L. Metody splayn-funktsiy [Methods of Spline Functions]. Moscow: Nauka Publ.,
1980. 352 p.
9. Fletcher K. Chislennye metody na osnove metoda Galerkina [Computational Galerkin Methods]. Translated from English. Moscow:
Mir Publ., 1988. 352 p.
10. Strang G., Strela V. Short wavelets and matrix dilation equations. IEEE Trans. Signal Processing. 1995, vol. 43, no. 1, pp. 108-115. DOI: 10.1109/78.365291
11. Stollnitz E.J., DeRose T.D., Salesin D.H. Wavelets for Computer Graphics: Theory and Applications. San Francisco: Morgan Kaufmann, 1996. 245 p.
12. Arandiga F., Baeza A., Donat R. Discrete multiresolution based on hermite interpolation: computing derivatives. Communications in Nonlinear Science and Numerical Simulation, 2004, no. 9, pp. 263-273. DOI: 10.1016/S1007-5704(03)00116-3
13. Ortega J.M. Introduction to Parallel and Vector Solutions of Linear Systems. New York: Premium Press, 1988. 305 p.
14. Starchenko A.V. (ed.) Vysokoproizvoditel'nye vychisleniya na klasterakh [High-performance calculations on clusters]. Tomsk: Tomsk State University Publ., 2008. 198 p.
