О феномене ложных корней систем линейных алгебраических уравнений при решении в действительных числах

Тезисы

1. При решении больших или плохо обусловленных систем уравнений часто возникает ситуация, когда прямая подстановка корней в уравнения дает прекрасную невязку для левых и правых частей, а на самом деле вместо ожидаемых корней получается чепуха. Покажем это на примере матрицы Гильберта 13*13, принцип формирования которой понятен из примера 3*3:

Хотя решение и очевидно (X₁=X₂=…=X_n=1), компьютер с ним не справится. Решая систему при двойной точности вычислений (Double, 15 десятичных разрядов), получим вместо единиц вектор неизвестных, показанный в таблице.
Применительно к задачам строительной механики корабля, матрица, похожая на матрицу Гильберта, получится, если первую силу приложить в центре пролета балки, вторую в центре между местом приложения первой и любой из опор, третью в центре между местом приложения второй и опорой и т.д. Кроме этого, под каждой из сил необходимо поместить упругую опору или перекрестную балку. Сама по себе задача расчета балки на регулярных упругих опорах является «плохой» в вычислительном плане из-за разностей близких величин. И она «ухудшается» многократно, если эти опоры нерегулярны и учащаются по мере приближения к крайней опоре. В то же время матрица Гильберта может быть задана абсолютно точно, поскольку выражается отношением целых чисел. Поэтому ее применяют для самого ответственного тестирования вычислительных программ.

Из таблицы видно, что строка с наиболее неудачным корнем дает нулевую невязку. Этот феномен хорошо известен специалистам, но плохо освещен в литературе. Нет прямых методов оценки погрешностей корней для диагностики "расхождения" решений при вполне удовлетворительных невязках. Все косвенные методы построены на анализе влияния "малых" возмущений правых частей или коэффициентов системы на изменение корней. Если "малые" возмущения приводят к "большим" изменениям корней, то ситуация признается "плохой".
Нет необходимости пояснять, что указанная диагностика напоминает настройку музыкального инструмента "на слух" мастера. А если вы еще не мастер, но уже понимаете, что прямая подстановка корней в уравнения не дает никаких гарантий их истинности? Тогда подойдет способ прямой проверки погрешностей корней, который так прост, что приходится удивляться, почему он не был открыт до сих пор. Этот метод получил название метода единичных корней.

2. Метод основывается на двух важных свойствах системы уравнений:
а) матрица коэффициентов имеет фундаментальные внутренние качества, не зависящие от значений правых частей, например, детерминант. Из этого следует, что если матрица "плохая", то она будет оставаться таковой для любых правых частей;
б) правые части, вычисленные для единичных корней, будут самыми точными из всех возможных ситуаций, поскольку умножение коэффициента системы на единицу не вносит никаких вычислительных погрешностей. Вычислив правые части, приступаем к решению и сравниваем потом полученные корни с единицей.
Хотя метод и прост, основывается он на неочевидных предположениях и требует серьезных размышлений. Именно простота метода рождает у читателя множество вопросов. Например, где гарантия, что он подойдет для конкретной задачи, где по физическому смыслу истинные корни должны быть разными? Но дело в том, что единичные корни не претендуют на роль «истинных» для физической задачи, их назначение – контроль правильности решения матрицы коэффициентов. Для этого необходимо идеальное соответствие правых и левых частей, достичь которого наиболее просто при единичных корнях. Более подробные доказательства приведены в основной работе [1].
Еще одно замечательное свойство единичных корней состоит в том, что они показывают нижнюю границу вычислительных погрешностей по каждому из корней. Это значит, что погрешности истинных корней для «физических» правых частей всегда будут хуже. Но насколько хуже? Показать теоретически то, что если единичный корень неверен в третьей значащей цифре, то примерно такая же погрешность и у «физического» корня, сложно. Зато очевидно, что если у единичного корня нет ни одной верной значащей цифры, то это же утверждение заведомо применимо к любому другому корню.

3. Что делать в том случае, если единичные корни показывают неудовлетворительные вычислительные погрешности? Ответ прост: независимо от того, насколько верны исходные данные, необходимо именно для них показать соответствующие корни, не выходя при этом за пределы установленных вычислительных погрешностей. Именно вычислительных, пока они не будут устранены, не может быть серьезных претензий ни к физической постановке задачи, ни к точности задания исходных данных. Это можно сделать, только повышая разрядность вычислений, другого пути нет. Как это осуществить практически? Есть специальные математические приложения, MathCAD, например. Есть специальные математические пакеты (профессиональные и нет) для языков программирования.
К сожалению, воспользоваться профессиональным пакетом автор не смог, так как не нашел бесплатного, с исходными кодами и с комментариями русском языке. Переходить в другую среду программирования (MathCAD) неудобно. Пришлось создать свой собственный пакет. Он уступает профессиональным по всем параметрам, зато бесплатный, с исходными кодами и с подробными комментариями. Написан он на языке Паскаль, его работа показана и проверена на примере решения системы уравнений.

4. Особенности обработки чисел современными процессорами таковы, что мантиссы чисел с повышенной разрядностью («длинных чисел») удобно представлять 32-х разрядными порциями. Результат расчета для 64-х битного числа (двойная точность, double) показан выше, среди корней есть настолько неудачные, что не содержат ни одной верной значащей цифры. Попробуем повторить расчет для числа из 96 бит. В этом случае получим единичные корни, верные до 10-ти первых значащих цифр.
В связи с этими результатами возникает вопрос, нет ли здесь расхождений с рекомендациями Лапласа, дополненных академиком А.Н. Крыловым [2]. Если говорить о точности корней с точки зрения вычислительного процесса, то нет. Что же касается физического содержания результата, то он, безусловно, не может быть «точнее» исходных данных. Как далеко можно заходить в длинных числах, предназначенных для расчета, при фиксированной длине исходных данных?
Решая систему для матрицы Гильберта 200*200, с точностью в 32*32 бит (приблизительно 308 десятичных разрядов) получаются верными 5 первых значащих цифр корней. Это притом, что длина чисел для расчета примерно в 16 раз превосходит длину исходных данных, которые для этого примера определены стандартом в 80 бит (extended). И это далеко не предел, если речь идет только о вычислительных погрешностях. Если задать точность в 32*33 бита, то верны 16 первых значащих цифр корней. Когда же наступает предел, ведь точность расчета нельзя повышать бесконечно при фиксированной точности исходных данных? Это можно установить только экспериментально, и одним из проявлений предела, может служить факт вырождения системы для высоких порядков матриц. Разумеется, что система не должна физически вырождаться ни при каких порядках матриц.

5. В основной работе приведена статья в полном объеме и три программы на Паскале (Delphi). В первой программе (приложение 1) осуществляется диагностика вычислительных погрешностей на примере матрицы Гильберта. Во второй и третьей реализован способ снижения погрешностей за счет привлечения арифметики длинных чисел. В приложении 2 реализована условно-десятичная арифметика, а в приложении 3 – двоичная. Это не одно и тоже, если речь идет о вычислениях с действительными числами, где постоянно приходится усекать мантиссу, округляя последнюю значащую цифру.

Полный объем работы с учетом трех программ в исходных кодах равен примерно 3500 строк. Поэтому здесь излагаются только основные положения, а электронную копию полного документа все заинтересованные читатели получат, обратившись к автору с письмом по адресу BLACK_EN@MAIL.RU.

( Основная работа, 53 страницы + программа в исходных кодах на Паскале (Delphi), 151кб )