Transcript В.Ю.Протасов (МГУ) Вычисление совместных спектральных характеристик линейных операторов A1, , Am - линейные операторы в Rd (1960) (1960) (1995) (1995) Joint spectral radius, Rota, Strang Lyapunov exponent, Furstenberg, Kesten,

В.Ю.Протасов (МГУ)
Вычисление совместных спектральных характеристик линейных операторов
A1,
, Am - линейные операторы в Rd
(1960)
(1960)
(1995)
(1995)
Joint spectral radius, Rota, Strang
Lyapunov exponent, Furstenberg, Kesten, (1968) Oseledec
p-radius, Lau, Wang
Lower spectral radius, Gurvits
The Joint spectral radius (JSR)
J.C.Rota, G.Strang (1960) -- Normed algebras
ˆ ( A1 ,
, Am )  lim
1/ k
max
k  d1 ,..., dk {1,..., m}
Ad 1
Ad k
Molchanov , Pyatnicky, Opoytsev, Micchelli, Prautzch, Dahmen,
Barabanov, V.Kozyakin,…(1988)
Levin, Dyn, …… (1989)
Linear switching systems
Subdivision algorithms
I.Daubechies, J.Lagarias (1991)
Wavelets
Совместный спектральный радиус (JSR)
, Am -- линейные операторы в Rd
A1,
ˆ ( A1 ,
, Am )  lim
1/ k
max
k  d1 ,..., dk {1,..., m}
Ad 1
Ad k
Геометрический смысл:
ˆ  1  существует норма
для которой
A2 M
M
A1M

в Rd g f
A i  1 при всех i  1, ... , m
Возьмём единичный шар в этой норме:
ˆ  1  существует центрально-симметричное выпуклое тело M  Rd ,g f
для которого Ai M  int M , i  1, ... , m
ˆ  inf {   0 |  1 A1 ,
,  1 Am  сжатия в некоторой норме}
Геометрический смысл JSR
x
Пусть m  2. Предположим, что семейство  A1, A2 неприводимо.
x  R dd  произвольная точка , x  0,
Ok ( x)  { Ad1 Adk x , d j  1, 2, j  1, , k } - орбита точки x.
 max
 |u|
| u  Ok 
JSR - показатель роста норм
A2 x
A1 x
O1 ( x) , O2 ( x) , O3 ( x) , ...
ˆ k
max Ad1 .... Adk
d1 ,..., dk

ˆ k
Приложения:
1960
1988-90
1991
1989-92
Rota, Strang (теория нормированных алгебр)
Барабанов, Козякин (динамическе системы с переключениями)
Daubechies, Lagarias, Cohen, Heil, …. (теория всплесков)
Micchelli, Prautzsch, Dyn, Dahmen, … (уточняющие схемы – теория приближений и
дизайн кривых и поверхностей)
Распределение случайных рядов (теория вероятностей),
Асимптотика бинарной функции разбиения Эйлера (комбинаторика, теория чисел),
Емкость кодов, оценка числа неперекрывающихся слов, теория графов, ....
Основные свойства
m  1. Для одного оператора ˆ (A) =  (A) = lim
k 
Ak
1/ k

max |  j |
j 1,..., d
Если операторы A1 ,... , Am коммутируют, либо их матрицы симметричны,
либо их матрицы -- верхне (нижне) треугольные, то
ˆ ( A1 ,..., Am )  max   ( A1 ) ,...,  ( Am ) 
В общем случае, однако ˆ ( A1 ,..., Am )  max   ( A1 ) ,...,  ( Am ) 
Всплески с компактным носителем (compactly-supported wavelets)
{  j ,k } j ,kZ полная ортонормированная система (ПОНС) в L2 ( R),
 j ,k ( x)  2 j / 2 (2 j x  k )
Пример 1. Система Хаара (Haar, 1909),  ( x) 
[0, 1/ 2]  [1/ 2, 1]

1, 1
А.Хаар (1909), В.А. Котельников (1933), К.Э.Шеннон (1949),.
1,0
1,1
1980-90: С.Малла, И.Мейер, И.Добеши, Ч.Чуи, А.Коэн, В.Дамен, и др.
I.Daubechies (1988) – всплески с компактным носителем.
Преимущества всплесков:
Локализация (компактные носители),
Быстрые алгоритмы вычисления коэффициентов,
Характеризация функциональных пространств
..
.
Обработка сигналов
Теория функций, теория приближений
f ( x) 
c
k , j Z
c j,k  ( f ,
Диф. Уравнения (Вейвлет-Галёркин метод, и др.)
j,k
j,k

) 

j,k
( x) ,
f ( x) 
j,k
( x) d x
Построение всплесков. Масштабирующие уравнения.
Для построения системы всплесков с компактным носителем нужно решить
масштабирующие уравнение (refinement equation)
– разностное функциональное уравнение со сжатием аргумента.
N
c
 ( x) 
k 0
c0 ,
, cN
k
 (2 x  k ) ,
- последовательность комплексных коэффициентов.
cN  (2 x  N )
( x)
c0  (2 x)
c1  (2x 1)
.......
Это – обычное разностное уравнение, но с двоичным сжатием аргумента.
Когда мастабирующая функция  найдена, всплеск-функция 
явно выписывается:
 ( x) 
N

k 0
(1)k ck 1  (2 x  k ),
где c0 , ... , cN - коффициенты масштабирующего уравнения.
Примеры систем всплесков
1. Всплески Хаара (1909)
Масштабирующее уравнение:
 ( x)   (2 x)   (2 x  1)
1
1

1
1
0
0
1
1
2. Всплески Шеннона-Котельникова
 ( x) 

sin  x
x
 ( x) 
(1933, 1949)
sin 2 x  sin  x
x
Носитель некомпактен!
3. Всплески Мейера (1986), Всплески Баттла-Лемарье (1987) Носитель некомпактен!
4. Всплески Добеши (1988)
 ( x) 
2
Второй всплеск Добеши. Масштабирующее уравнение:
1 3
3 3
3 3
1 3
 (2 x) 
 (2 x  1) 
 (2 x  2) 
 (2 x  3)
4
4
4
4
2

0
3

3
0
Что известно о масштабирующих уравнениях ?
 ( x) 
N
c
k 0
k
 (2 x  k ) ,
Если есть решение с компактным носителем, и
N
Обратно, если
c
k 0
k
  ( x) d x
 0
N
то
c
k 0
k
 2
 2 , то есть решение с компактным носителем. Оно единственно
с точностью о домножения на константу и сосредоточено на отрезке [0, N].
( x)
Но только в обобщённых функциях из
S '(R )
0
ˆ ( )

m( / 2) ˆ ( / 2) ,
ˆ ( )


N
1 N
m( ) 
c k e  2  i k

2 k 0
m ( 2
j
)
j 1
Масштабирующая функция не бывает бесконечно-гладкой
  C N (R)
Примеры масштабирующих уравнений
Примеры 1. c0  c1  1
Тривиально:
 ( x)   (2 x)   (2 x 1)
Пример 2.
Пример 3.
c0 
c0 


0
1
1
, c1  1 , c 2 
2
2


1
3
3
1
, c1  , c 2  , c 3 
4
4
4
4


0
1
2
0
3
Решение неустойчиво !
Малые изменения коэффициентов могут приводить к резким изменениям решения:
Пример 4.
c0 
1

2
1

Tо же с примером c 0 
4
c1  1  
3
c1   
4
c2 
1
2
c2 
3
4
 чисто сингулярно.
c3 
1
4
 чисто сингулярно.
“ Типичная ” масштабирующая функция и всплеск-функция
Пример 5
 ( x)

1
2
1
2
 (2 x)   (2 x  1)   (2 x  2)   (2 x  3)
3
3
3
3
 ( x)  log 2 3  1.58...
(максимальная гладкость)
 ( x)
 ( x)  log 2 1.5  0.58...
(минимальная гладкость)
0

 sup    0 |
(изломы во всех двоичнорациональных точках)
3
|  ( x  h)   ( x) |  C h  for all x, h 
показатель гладкости (показатель Гельдера)
  log 2 1.5  0.58...


Непрерывна, но не дифференцируема
Тем не менее, она дифференцируема почти всюду
 ( x)  sup    0 |
|  ( x  h)   ( x ) |  C h 

Локальная гладкость в точке x
Почти во всех точках  ( x)  log 2 2.25  1.17...
Следовательно,  '( x)  0 п.в.
Фрактальная природа всплесков. Переменная локальная гладкость.
Как определить, будет ли масштабирующая функция непрерывной ?
I.Daubechies, D.Lagarias, 1991
A.Cavaretta, W.Dahmen, C.Micchelli, 1991
 C (R)  ˆ (T0 |W , T1 |W )  1
Более того,    log 2 ˆ (T0 |W , T1 |W )
C.Heil, D.Strang, 1994
ˆ ( A0 , A1 )  lim
k 
1/ k
max
d1 , , dk
Ad 1
Ad k
Пример.
N  4, c0 , c1 , c2 , c3 , c4
 c0

c2

T0 
 c4
T0 , T1 - N  N  матрицы (2-блочные тёплицевы матрицы),

0
(Ti ) j k  c2 j  k 1i
 c1

c
T1   3
0

0
0
0
c1
c0
c3
c2
0
c4
c0
0
c2
c4
c1
c3
0
0
0

0
c1 

c3 
0

c0 
c2 

c4 
Как вычислить или оценить JSR ?
Перебором (  по определению) Daubechies, Lagarias, Heil, Strang,
Эспоненциальная сложность. Heil, Strang (1994) для вычисления JSR специальных 2  2-матриц
с относительной погрешностью  = 0.05 перебрали все произведения до длины k  19.
G.Gripenberg (1996) - ``branch and bound'' алгоритм.
Разумный перебор. Часто очень эфективен. Но теоретически плох.
Сходимость к величине JSR при растущем к
очень медленная.
Причина медленной сходимости:
1/ k
max
d1 ,..., d k {1,..., m}
Ad 1
Ad k

ˆ ,
k 
Выбранная норма в R d d может быть слишком далека
от той, в которой все операторы - сжимающие.
C
, где константа C - отношение двух норм.
k
Она может быть велика.
Скорость сходимости
Инвариантные нормы
Теорема 1 (Н. Барабанов, 1988)
a) Для неприводимого семейства операторов A1 ,..., Am существует   0
и норма

max
, для которой при любом x  R d

A 1 x , ... , A m x

=  x
b) Для любой такой нормы имеем   ˆ ( A1 ,..., Am ).
Независимо был установлен ‘’двойственный’’ факт:
A2M
Теорема 2 (A.Дранишников, С.Конягин, В.Протасов, 1996)
M
a) Для неприводимого семейства операторов A1 ,..., Am существует   0
A1M
и симметричное выпуклое тело M (инвариантное тело) такое, что
def
A M  Conv ( A1 M ,..., Am M )   M
b) Для любого такого тела   ˆ ( A1 ,..., Am ).
def
A M  Conv ( A1 M ,..., Am M )
Двойственность: M = B* (поляра к B), где B - единичный шар в Барабановской норме
для сопряженных операторов A1* ,..., Am* (F.Wirth, E.Plishke, 2005)
Все существующие методы вычисления JSR основаны на
приближении инвариантного тела
(многогранниками, эллипсоидами, и т.д.)
Branch-and-bound method («разумный перебор»)
Grippenberg (1996), Hechler, Mossner, Reif (2009)
Ellipsoidal norms,
Daubechies, Lagarias (1991),
Ando, Shih (1998), Blondel, Nesterov, Theys (2004)
Kronecker lifting method (приближение инвариантной нормы тензорными
произведениями) , Protasov (1997), Zhow (1998), Blondel, Nesterov (2005)
Sum-of-squares method (приближение инвариантной нормы суммами квадратов
полиномов), Parrilo, Jadbabai (2008)
Conic programming approach, Protasov, Jungers, Blondel (2010)
Все методы – для невысоких размерностей (до 5, реже – до 10)
Точность невысока. Относительная погрешность – от 2 до 20 %
Отрицательные результаты о сложности задачи вычисления JSR:
Blondel, Tsitsiclis (1997-2000).
Задача приближённого вычисления JSR для рациональных матриц NP-сложна.
Задача определения: верно ли, что JSR строго меньше 1
(для рациональных неотрицательных матриц) алгоритмически неразрешима,
начиная с размерности d = 47.
1
Не существует алгоритма, полиномиального по размерности d и по точности

для приближения JSR с относительной погрешностью

Sometimes easier to prove more
George Polya «Mathematics and Plausible Reasoning» (1954)
When trying to prove something, often a good strategy is to try to prove more.
(Если не получается что-то доказать, можно попробовать доказать больше).
Если не получается хорошо приблизить JSR, можно попробовать …
вычислить его точно.
Метод точного вычисления совместного спектрального радиуса
(N. Guglielmi, V.Protasov, 2013).
Применим к любым размерностям (эффективен на практике в
размерностях до 20, для неотрицательных матриц – до 100)
Сводит задачу вычисления JSR к поиску корня некоторого полинома
Позволил решить ряд открытых проблем теории чисел и комбинаторики
Понятие экстремальной нормы
Определение. Норма
.
является экстремальной, если
Aj
 ˆ , j  1,..., m .
Для экстремальной нормы сходимость осуществляется в один шаг:
1/ k
max
d1 ,..., dk {1,..., m}
Ad 1
Ad k
 ˆ ,
k 
Если M - единичный шар экстремальной нормы, то A j M  ˆ M ,
A2 M
A1M
M
j  1,..., m.
Будем строить единичный шар экстремальной нормы в качестве многогранника M .
Оказывается, что такая норма существует для большинства семейств матриц.
Наблюдение 1. Для приводимого семейства задача вычисления JSR сводится к
Нескольким аналогичным задачам в меньших размерностях. Таким образом,
предполагаем, что семейство неприводимо.
Наблюдение 2. Если произведение П максимальное, то его максимальный
собственный вектор v должен быть крайней точкой множества M.
Если M – многогранник, то v -- его вершина.
Итак, максимальные собственные векторы произведения П и всех его циклических
перестановок -- вершины M.
Qv
Наблюдение 3. Критерий остановки:
v
Лемма. Пусть  ( ) = 1 и v* - максимальный
собственный вектор * . Если существует другое
произведение Q, для которого (v*, Qv) > (v*, v),
то 
-
не
максимальное.
v
Qv
Алгоритм точного вычисления JSR (N.Guglielmi, V.Protasov, 2011)
Берем максимальный собственный вектор v1 of   Ad1
Полагаем v j  Ad k  j2
j  2,
Adk v1 ,
vk
, k.
…..
Ad2
Ad1
v1
Adk .
Adk
Adk2
v2
v3
Adk1
vs  A j v1
v p  Ai v q
‘’Мертвые’’ ветви
Каждый раз проверяем, будет ли новая точка принадлежать выпуклой оболочке
предыдущих точек (ЛП задача).
Алгоритм завершается, когда не появилось ни одной новой вершины.
Инвариантный многогранник M – выпуклая оболочка всех точек, построенных алгоритмом
Пример 1. Задача о плотности единиц в ромбе Паскаля:
(S.Finch, P.Sebah, and Z.-Q.Bai, 2008)
Плотность единиц в ромбе Паскаля порядка n -- между C nlog2  ( A1 , A2 ) и
C nlog2 ˆ ( A1 , A2 ) , где
Известно, что ˆ ( A1, A2 ) = 2.
Относительно  ( A1 , A2 ), выдвинута гипотеза, что он равен
1+ 5
= 1.6180...
2
На самом деле,
 ( A1 , A2 )    A A2
3
1
(алгоритм работает несколько секунд)

3 1/ 6
 1.6376...
Выбираем   A13 A23
Инвариантный многогранник M1 имеет 8 вершин.
Пример 2. Асимптотика бинарной функции разбиения Эйлера.
Бинарная функция разбиения Эйлера bd (k ) -- это количество разложений
k  d0  21 d1  22 d2 
2m1 dm1 ,
где d j {0,1,
, d 1}
Как мы знаем, b2 (k )  1. Для d  3 нужно оценить рост bd (k ) при k  .
b2 (k )
 1
(L.Euler, 1728)
b3 (k )  s(k  1) (Stern, 1858)
b4 (k ) 
k / 2  1
(Klosinsky, Alexanderson, Hillman, 1984)
Какова асимптотика величины bd (k ) при k   ?
L.Euler (1728), A.Tanturri (1918), K.Mahler (1940), N.de Bruijn (1948)
L.Carlitz (1965), D.Knuth (1966), R.Churchhouse (1969), B.Reznick (1990)
Ответ:
b2 r (k )  C k  , где   log 2 r
(B.Reznick, 1990)
C k 1 
b2 r 1 (k )  C k 2 , где
1  log 2 ( A1 , A2 ),  2  log 2 ˆ ( A1 , A2 ) (V.Protasov, 2000)
A1 , A2 это (d  1)  ( d  1)  матрицы из нулей и единиц:
 1, если 2  2k  j  i  d  1
( Ai ) j k  
 0, иначе.
1

0

A1 
0

0
1
1
1
0
Пример. При d = 5 :
1
1
1
1
0

0
1

1
1

1
A2  
0

0
Алгоритм вычисляет точные значения для d  100.
Оказывается, что для всех d имеем
либо ˆ 
( A 1 A2 ) ,   ( A1 ) , либо:  
( A 1 A2 ) , ˆ  ( A1 )
Для размерности 50 программа работает 5 минут,
для размерности 100 -- около 20 минут
1
1
1
1
0
1
1
1
0

0
0

1
Функция разбиения Эйлера для троичного разложения:
Выбираем   A1 A3
Экстремальный многогранник M3 имеет 16 вершин.
Пример 3. Асимптотика числа слов двоичного алфавита без перекрытий.
Задача сводится к вычислению JSR и LSR двух 20x20-матриц.
  2 (Blerk, 1988) ,   2.226 (Kobayashi, 1988)
ˆ  2.584 Этот результат последовательно улучшался: 2.226    ˆ  2.584
Kforu (1988), Kobayashi (1988), Cassaigne (1993), Lepisto (1995)
    ( A1 A210 ) 
1/11
 2.41756...
;
ˆ    ( A1 A2 )
Программа работает 8 минут
1/ 2
 2.51793...
Вычисление JSR для случайных пар матриц
Вычисление JSR и LSR для случайных пар булевских матриц размерности d = 100.
Условия конечной сходимости алгоритма
Определение. Произведение   Adk
Ad1 называется доминирующим, если  ( )=1,
и существует q  1 такое, что  () < q для всех остальных произведений   Adn
не являющихся степенями , или степенями его циклических перестановок.
доминирующее
максимальное
Теорема 1. Алгоритм сходится за конечное время тогда и только тогда
когда произведение  доминирующее.
Ad1 ,
Спасибо!