Унитарная,минимально удаленная от нескольких других унитарны

[zgemm]

Заголовок: Найти такую унитарную, которая была бы минимально удалена от нескольких других унитарных

Добрый день и с наступившим Вас новым годом, пусть этот квадратичный год возведет все Ваши загаданные желания в квадрат и они сбудутся!

Задачка, что я к Вам обращаюсь - в заголовке, то есть у меня есть набор унитарных матриц $Q_n \in \R^{m \times m}, n=1,...,N$, и я ищу такую унитарную $P$ которая бы была бы минимально удаленной (наверное по Евклидовой норме) от всех этих.

То есть примерно так:

$$\min_P \sum_{n=1}^N ||Q_n - P||_2^2 + \lambda ||P^* P - I||_2^2$$

где $\lambda$ - множитель Лагранжа.

Если попробовать явно записать частные производные и приравнять их к нулю, то получится полиномиальная задача с лагранжевым множителем и очень не очевидно какая тут будет сходимость у методов, типа квазиньютона например.

Наверное можно еще записать

$$E(\lambda) = \min_P \sum_{n=1}^N ||Q_n - P||_2^2 + \lambda ||P^* P - I||_2^2$$

и стартовать такую минимизацию для фиксированного и довольно маленького $\lambda$, а потом его увеличивать и увеличивать, пока не состоится унитарность $P$ с любой на перед заданной точностью.

Других решений в голову не приходит.

Вдруг у Вас будут еще какие-то идеи, пожалуйста, посоветуйте!

Спасибо!

[Евгений Машеров]

Идея из разряда дурацких:

Среднее геометрическое унитарных матриц. В смысле перемножить и вычислить степень $\frac 1 n$
При перемножении унитарность сохраняется, а корень n-ной степени можно попробовать через матричный логарифм и экспоненту.

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

[dgwuqtj]

А ещё можно взять среднее арифметическое и ортонормировать её столбцы. Это как первое приближение может сгодиться.

[Евгений Машеров]

Ну, вообще для унитарных матриц есть представление $U=e^{iH}$, где H - эрмитова матрица.
В Матлабе есть функция logm и expm. Возможно, так получится - найти $iH_j=\log U_j$, усреднить $\bar{iH}=\frac 1 n \sum iH_j$ и затем получить $P=e^{\bar{iH}}$

[zgemm]

Спасибо Евгений Машеров! Про среднее геометрическое - не подумал, идея довольно красивая, надо попытаться ее далее развить.

Спасибо dgwuqtj! Усреднить и потом ортонормировать уже пробовал. Если ортонормировать Граммом-Шмидтом, часто получается лажа. Грубо говоря пусть у меня есть матрицы, они в моей задаче обычно довольно похожи друг на друга. Когда беру среднее - да, очевидно получается не унитарная. Если ортонормировать столбец за столбец, то первые столбцы меняются мало, а вот те, что в конце - меняются уже довольно сильно. В результате, подставляя такую матрицу в исходное выражение для невязки у меня невязка будет больше, чем если я возьму любую из исходных унитарных.

Заметил, что

$$\min_p \sum_n ||Q_n - P||_2^2 + \lambda ||P^* P - I||_2^2 =$$

$$\min_p \sum_n ||Q_n||_2^2 - 2 \sum_n \operatorname{tr}(Q_n^*P) + \sum_n ||P||_2^2 + \lambda ||P^* P - I||_2^2 =$$

$$ \sum_n ||Q_n||_2^2 + N \min_p || \frac{1}{N} \sum_n Q_n - P||_2^2 + \frac{\lambda}{N} ||P^* P - I||_2^2 =$$

то есть на самом деле задача сводится к поиску унитарной, наиболее близкой от некоторой на перед заданной (то есть средней арифметической).

[Евгений Машеров]

А для ортогонализации среднего не сгодится ли полярное разложение?

[Евгений Машеров]

Примерно так. Наша "усреднённая" матрица A может быть представлена в виде $ A=SU$, где S - эрмитова, и U - унитарная. Если исходная матрица уже унитарная, то $ S=I $, а если нет, то матрица S не единичная. Она "выбирает" неунитарность матрицы A, и если взять только U - возможно, это ответ.
В принципе, можно оба подхода испытать - со среднегеометрическим (или, что то же, с усреднением логарифмов матриц) унитарность сразу обеспечивается, но средний квадрат отклонения может быть выше, а если начинать со среднего арифметического, то приближение самое лучшее, но при "унитаризации" могут быть потери и в целом хуже.

[zgemm]

Не, там все будет проще. Делаем для этой матрицы сингулярное разложение $S = UDV^*$. Очевидно, что теперь все можно заменить так, что вместо $S$ мы работаем с диагональной $D$. Так как среднее матриц не может увеличить сингулярное значение, в этой диагональной все элементы на диагонали не больше единицы. Следовательно, решение надо искать только на классе диагональных, а в этом классе есть только единичная. Тогда заменяем $D$ на единичную и получаем искомую в виде $P=UV^*$.

PS: почему мне это сразу в голову не пришло - хз, наверное праздники и сопутствующие мероприятия сказываются. Всем с наступившим!!!

[Евгений Машеров]

Ну, собственно, полярное разложение это та же идея.
$UDV^*=UV^*VDV^*=(UV^*)(VDV^*)=WS$
W унитарная
S эрмитова
Но интересно было бы и со средним геометрическим попробовать, через усреднение матричных логарифмов.

[zgemm]

Евгений Машеров в сообщении #1668377 писал(а):

Ну, собственно, полярное разложение это та же идея.

так правильно, я из-за Вашего сообщения в эту сторону посмотрел, спасибо большое!!! Другое дело, что с сингулярным разложением все доказуемо и получается аналитически, то есть пока мы из $S$ в $D$ преобразуем, все нормы сохраняются, а далее, замечаем, что любое сингулярное у $S$ не может быть больше, чем сингулярное от любой исходной, что приводит к тому, что на диагонали $D$ стоят не отрицательные числа не более единицы, а для такой матрицы ближайшая по евклидовой норме унитарная будет единичной.

[Geen]

zgemm в сообщении #1668415 писал(а):

стоят не отрицательные числа не более единицы

а это важно?

[zgemm]

Geen в сообщении #1668423 писал(а):

а это важно?

похоже, Вы правы, и это не важно, спасибо! Я почему-то вчера думал, что если взять два на два матрицу с одним диагональным элементом больше единицы и одним - меньше, то для такой матрицы ближайшая по евклидовой норме будет не диагональная, а сейчас перепроверил и понял, что ошибался.

PS: да праздники все-таки плохо действуют на работу...

[Евгений Машеров]

Интересно было бы узнать, получилось ли.

[zgemm]

Евгений Машеров в сообщении #1669184 писал(а):

Интересно было бы узнать, получилось ли.

вроде же через среднее и SVD доказалось (смотрите обсуждения выше), а это заимплементировать - это же всего несколько строк на С/С++ с лапаком и, очевидно, я это как только доказал, так сразу и запрограммировал, и уже несколько дней как использую.