Норма невязки

[Andrey_Kireew]

При итеративном уточнении обратной матрицы методом Шульца

$$B^{-1}_{k+1}=B^{-1}_{k}(2E-B\cdot B^{-1}_{k})$$

рано или поздно, а на практике - уже через 2-3 итерации, достигается предельная точность, с смысле нормы невязки

$$\delta =||E-B\cdot B^{-1}_{k}||$$

Так как объективно, в исходных данных для расчёта имеется лишь матрица $B$ и eps, разумно предположить, что именно они и определяют эту остаточную невязку. Причём, у меня возникает предположение, что

$$\delta \to ||B||\cdot eps $$

Верно ли это предположение?

[Andrey_Kireew]

В общем, возникли такие соображения:

$$(B+\Delta_B)(B^{-1}_k+\Delta_{IB})=E+\Delta$$

Раскрываем скобки

$$ B\cdot B^{-1}_k+\Delta_B \cdot B^{-1}_k+B\cdot\Delta_{IB}+\Delta_{IB}\Delta_{B}=E+\Delta$$

,
упрощаем, учитывая, что произведение $\Delta_{IB}\Delta_{B}$ пренебрежимо мало

$$ \Delta=\Delta_B \cdot B^{-1}_k+B\cdot\Delta_{IB}$$

.

Предположим, что

$$||\Delta_B ||\leqslant ||B||\cdot eps$$

и

$$||\Delta_{IB} ||\leqslant ||IB||\cdot eps$$

тогда, получается

$$ ||\Delta||\leqslant 2\cdot cond(B) \cdot eps$$

Для проверки, я вычислял следующее отношение

$$\frac{||\Delta||}{2\cdot cond(B)} $$

для разных матриц после установления итеративного процесса.
В результате, получились значения, лежащие в диапазоне $1,555\cdot 10^{-17}...3,1920\cdot 10^{-17}$, т.е. максимально различающиеся всего в 2 раза.
Смущает только одно. Эти значения примерно в 8 раз меньше чем eps для чисел double, которое, если мне не изменяет память равно $eps=2,22\cdot 10^{-16}$.

Возникает интуитивное чувство, что в вышеприведённой формуле, вместо коэффициента 2 должен быть коэффициент 1/4.

Нет ли в моих рассуждениях ошибок?

[vpb]

Andrey_Kireew

Бывают три разных оценки для ошибки: (а) теоретически доказуемая оценка, (б) результат какой-то
прикидки на пальцах, неизвестно из каких соображений, и (в) результат численного эксперимента. И они совпадать совсем не обязаны. Это не только для невязки при вычислении обратной матрицы по Шульцу, а и во всяких других ситуациях.

Теоретически, обратную матрицу можно вычислить с невязкой порядка $Cn^2\operatorname{cond}(A) \varepsilon$ (емнип). А практически всегда получается лучше (например из-за того, что в ходе вычисления возможные ошибки округления между собой сокращаются случайным образом). От результатов численного
эксперимента, как видите, "дистанция огромного размера" (в $\sim n^2$ раз).

Я когда-то тоже размышлял, как можно поточнее вычислять обратную матрицу. Не помню уж, чем там
размышления кончились. Главная философия насчет обратной матрицы, согласно книжкам, такая: почти всё, для чего её можно использовать, быстрее и проще делается без вычисления обратной матрицы. Например, когда решаем систему $Ax=b$, там же обратную матрицу $A^{-1}$ не вычисляем, не так ли ?

Про "остаточную невязку", строго говоря, говорить нельзя. Разве что про оценку для нее. Потому что
сначала результат итерации Шульца приближается к обратной матрице, потом попадает в некоторую ея
окрестность, и потом уже в этой окрестности топчется до бесконечности туда-сюда-и-в-стороны.
(Сам я никогда экспериментов такого рода так не делал (оно мне и не надо), но из общих соображений
должно быть именно так.) Если бы все вычисления делались точно, то, понятно, точность бы нарастала. Но в реальных вычислениях с плавающей точкой такого, ясно, быть не может.

[Andrey_Kireew]

vpb в сообщении #1452931 писал(а):

обратную матрицу можно вычислить с невязкой порядка $Cn^2\operatorname{cond}(A) \varepsilon$ (емнип)

уточните пожалуйста, что означает в этой формуле n? размер матрицы?

-- 09.04.2020, 00:31 --

vpb в сообщении #1452931 писал(а):

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

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

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

-- 09.04.2020, 00:35 --

vpb в сообщении #1452931 писал(а):

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

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

[vpb]

Andrey_Kireew в сообщении #1452938 писал(а):

размер матрицы?

Именно так.

Andrey_Kireew в сообщении #1452938 писал(а):

Вот этот предел мне оценить и нужно.

Это я, увы, не могу сказать, как сделать. Возможно, есть работы, в которых даются соответствующие оценки (типа " при любых начальных условиях начиная с некоторого момента невязка не превосходит $n\operatorname{cond}(A)\varepsilon$". Тут я фантазирую, как на самом деле, не знаю (и есть ли такие работы вообще). ) Могу подсказать, как их можно поискать. Найдите хоть какую статью про итерацию Шульца, из нее извлеките ссылку на исходную работу Шульца, потом вбейте её в гугл сколар и посмотрите, в каких статьях на нее ссылаются. Или решайте, когда оборвать итерации, чисто интуитивно, по тому, что численные эксперименты подсказывают.

[Andrey_Kireew]

vpb в сообщении #1452951 писал(а):

Именно так

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

vpb а что Вы можете сказать на счёт вышеприведённых формул? На Ваш взгляд они справедливы, или ошибочны?
кстати, может быть я неудачно выбрал нормы матриц? Сейчас у меня вычисляются максимальные суммы модулей в столбцах, просто потому, что их вычислить проще всего

[vpb]

Andrey_Kireew в сообщении #1452956 писал(а):

что то мне подсказывает, размер в этой формуле фигурировать никак не должен. Для прямых методов - возможно, но не для итеративных.

Этого не может быть. Дело в том, что уже при умножении матрицы на вектор в оценке ошибки фигурирует размер.

Andrey_Kireew в сообщении #1452956 писал(а):

а что Вы можете сказать на счёт вышеприведённых формул? На Ваш взгляд они справедливы, или ошибочны?

Вы имеете в виду формулу $\|\Delta\|\leq 2\varepsilon \operatorname{cond}(B)$, вероятно. По результатам численного эксперимента --- вероятно, справедлива. А в строгом смысле скорее всего нет.

Andrey_Kireew в сообщении #1452956 писал(а):

кстати, может быть я неудачно выбрал нормы матриц? Сейчас у меня вычисляются максимальные суммы модулей в столбцах, просто потому, что их вычислить проще всего

Я имел в виду спектральную норму ... или как она называется? Короче, индуцированную обычной евклидовой нормой на векторах. Та, о которой Вы писали ---тоже индуцированная, но другой нормой на векторах (то ли 1-нормой, то ли $\infty$-нормой, сейчас не соображу).

[Andrey_Kireew]

vpb в сообщении #1452965 писал(а):

По результатам численного эксперимента --- вероятно, справедлива. А в строгом смысле скорее всего нет

Так весь и вопрос то в этом. В каком месте моих рассуждений закралась ошибка? я и сам подозреваю, что где то она есть ...

А на счёт размеров, я категорически несогласен. Число обусловленности с увеличением размеров, да - растёт (где то даже встречал формулу), но оно уже заложено в формуле. Или Вы сомневаетесь в том, что норма ошибки округление матрицы?

$$||\Delta_B||\leqslant ||B||\cdot eps $$

ведь размеры могут всплыть только в этом месте.

Признаться, я и сам сомневаюсь в этой формуле, но пока не знаю в чём загвоздка.

[Andrey_Kireew]

Проверил всё ещё раз, оказалось, что формула

$$||\Delta||\leqslant 2\cdot cond(B)\cdot eps$$

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

По поводу размеров - я тоже разобрался. Действительно, размеры матрицы входят в оценки точности решения СЛАУ прямыми методами, причём, в зависимости от конкретного метода, входят по разному. В случае итерационных методов, размеры матрицы, как я и предполагал, непосредственно не влияют на предельную точность решения, только - посредством возрастания числа обусловленности. Именно поэтому, в случае больших матриц прямые методы сильно уступают итерационным по точности, и практически не применяются.

[vpb]

Пусть $B$ --- невырожденная матрица, элементы которой --- числа с плавающей точкой, $B^{-1}$ --- обратная матрица (точная), и $C$ --- матрица, полученная из $B^{-1}$ поэлементным округлением. Тогда, несомненно, $\|C-B^{-1}\|\leq\varepsilon\|B^{-1}\|$. Отсюда $\|BC-E\|\leq\varepsilon\|B\|\|B^{-1}\|=\varepsilon\operatorname{cond}(B)$. Но, в результате процесса Шульца матрица $C$ отнюдь не получается. И вообще, какова оценка для $\|\widehat C-B^{-1}\|$, где $\widehat C$ --- матрица, получающаяся после более-менее долгого применения процесса Шульца, непонятно.

-- 13.04.2020, 21:43 --

vpb в сообщении #1454302 писал(а):

долгого применения процесса Шульца, непонятно.

Уточню. Я когда-то размышлял, у меня получилось в уме и на пальцах оценка типа $n^2\varepsilon\|B^{-1}\|$, не помню подробности. Что в литературе по этому поводу написано, не знаю. Т.е. может быть в науке оно и известно.

[Andrey_Kireew]

vpb в сообщении #1454302 писал(а):

Но, в результате процесса Шульца матрица $C$ отнюдь не получается

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

В общем, действительно, эта формула оценивает минимальную остаточную невязку, которая будет присутствовать в самом худшем случае, из-за ошибок округления. Будет ли она гарантированно достигнута в результате процесса Шульца, собственно и есть главный вопрос. Вернее - вопрос в том, какая невязка будет гарантированно достигнута. Лично я не вижу причин, почему не будет. Хотя это лишь интуитивное предположение, подкреплённое численным экспериментом. Но и Ваше утверждение точно такое же голословное.
А формула с квадратом размеров меня просто удручает. Я проводил вычисления с разными матрицами, в основном большими - до нескольких тысяч. Да, предельная невязка, оцениваемая по моей формуле получается в 8 раз (а если убрать двойку, то в 4 раза) больше, того что наблюдается. Какой либо заметной зависимости от размеров я не обнаружил.
Но оценки, полученные по Вашей формуле отличаются от действительности на много порядков. Спрашивается, какой толк от таких оценок? У меня такое подозрение, что они просто напросто неверны.

[vpb]

Да, действительно, разрыв между теоретически доказуемой оценкой и той, которая обнаруживается экспериментально, довольно большой. Но в вычислительной математике это дело обычное, се ля ви. Я когда-то потратил некоторое (небольшое, впрочем) количество усилий, думая, нельзя ли получить лучшую оценку, чем с $n^2$, типа того. Не преуспел, и желания возвращаться к этим размышлениям как-то нет.

Andrey_Kireew в сообщении #1456079 писал(а):

У меня такое подозрение, что они просто напросто неверны.

Неверны они быть не могут, так как они доказаны. А того, что оне неулучшаемы, я не утверждал.

Могу рекомендовать кой-какие источники. Именно, про численное обращение матриц кое-что написано в таких двух фундаментальных книгах, как Уилкинсон, Алгебраическая проблема собственных значений, и Higham, Accuracy and stability of numerical algorithms. Во второй целая глава посвящена.

[Andrey_Kireew]

Видимо я не корректно выразился, действительно, неверны они быть не могут, ведь они в $n^2$ больше. Вернее сказать - они недостаточно хороши.

Тем не менее, в книге Амосов А.А и др. Вычислительные методы для инженеров приводится оценка погрешности прямых методов решения СЛАУ

$$\delta < f(n)cond(B) eps$$

В случае QR разложения f(n)=6n, для метода Гауса в выбором ведущего элемента - теоретически больше, но сказано, что еще не обнаружено такой матрицы, для которой f(n)>n. Так, что погрешность растёт не с квадратом размера, а линейно. И это для прямых методов. Понятное дело, что итеративные методы всегда точнее. Поэтому я только укрепляюсь в своей уверенности, что погрешность у них от размера матрицы практически не зависит.

Уточните пожалуйста, Ваша оценка действительно получена применительно к итеративным методам?