В. В. ВОЕВОДИН ВЫ ЧИСЛИТЕЛЬНЫЕ ОСНОВЫ . ЛИНЕЙНОЙ АЛГЕБРЫ Допущено Министерством высшего и среднего специального образования СССР в качестве учебного пособия для студентов вузов, обучающихся по специальности «Прикладная математика» ИЗДАТЕЛЬСТВО «НАУКА» ГЛАВНАЯ РЕДАКЦИЯ ФИЗИКО-МАТЕМАТИЧЕСКОЙ ЛИТЕРАТУРЫ Москве 1977 518 В 63 УДК 519.95 Вычислительные основы линейной алгебры. В. В. Воеводин. Главная редакция физико-матема- тической литературы изд-ва «Наука», М., 1977. В книге последовательно изучаются ошибки округ- ления элементарных арифметических операций, их про- исхождение, свойства и влияние, на вычислительные процессы, впервые рассматриваются вероятностные свойства ошибок округления. Описываются основные . численные методы, связанные с решением систем, вы- числением определителей, решением полной и частичной проблем собственных значении. Особое внимание уде- ляется изучению устойчивости и оценкам влияния ошибок округления. Книга рассчитана на студентов вузов, изучающих прикладную математику, и будет полезна всем лицам, решающим задачи алгебры на ЭВМ. Библиографических названий 9. 20204—122, ^ © Главная редакция В S-11 физико-математической литературы им\и^)-11 издательства «Наука», 1977 ОГЛАВЛЕНИЕ Предисловие ............................. 5 Глава I Математические особенности машинной арифметики § 1. Позиционные системы счисления .............. 11 § 2, Округление чисел ....................... 15 § 3. Фиксированная и плавающая запятая ........... 20 § 4. Особенности представления чисел на ЭВМ ........ 23 § 5. Арифметические операции .................. 27 § 6. Порядок выполнения операций ............... 32 § 7. Запись на машинно-независимых языках ......... 36 § 8. Суммарный эффект влияния ошибок округления ..... 42 Глава II Теория возмущений в линейной алгебре § 9. Сведение к простым матрицам ............... 46 § 10. Невырожденные матрицы .................. 52 § 11. Непрерывность корней алгебраического многочлена ... 55 § 12. Локализаций собственных значений ............ 61 § 13. Клеточно-диагональные матрицы .............. 66 § 14. Матрицы общей структуры ................. 71 § 15. Сингулярное разложение .................. 73 § 16. Проекции псевдорешения .................. 78 § 17. Нормальное псевдорешение ................. 84 Глава III Вспомогательные алгебраические операции § 18. Преобразование вращения ................>•. 89 § 19. Последовательность преобразований вращения ...... 96 § 20. Преобразование отражения ................. 103 § 21. Последовательность преобразований отражения ...... Ill . § 22. Сравнение точности преобразований вращения и отражения 114 § 23. Двухсторонние унитарные преобразования ........ 117 § 24. Неунитарные преобразования ................ 122 § 25. Ортогонализация . . . . . . . . . . . . . . . . . . . . . . . 128 Г 4 ОГЛАВЛЕНИИ Глава IV Прямое разложение матрицы на множители § 26. Матрицы специального вида ................ 137 § 27. Теоретические основы разложения ............. 141 ' § 28. Разложение на треугольные множители .......... 146 § 29. Компактная схема ...................... 155 § 30. Разложение на унитарный и треугольный множители ... 161 § 31. Разложение прямоугольных матриц ............ 164 § 32. Унитарно подобное разложение .............. 167 ' § 33. Некоторые замечания .................... 172 § 34. Сравнительная характеристика разложений ......... 175 Г лава V Решение систем линейных алгебраических уравнений § 35. Системы специального вида ................. 179 § 36. Решение систем с невырожденными матрицами ...... 184 § 37. Системы с матрицами полного ранга ............ 190 . § 38. Уточнение решения ..................... 197 § 39. Особенности решения неустойчивых систем ........ 205 § 40. Системы с двухдиагональными матрицами ......... 211 § 41. Тактика решения систем общего вида ........... 215 § 42. Некоторые замечания .................... 223 Глава V/ Решение проблемы собственных значений $ 43. Метод вращении ....................... 229 S 44. Метод бисекций ....................... 237 § 45. Q^-алгорифм ......................... 249 § 46. Ускорение QR-алгорифма .................. 255 ^ 47. Определение собственных векторов ............. 2G4 § 48. Особенности вычислений .................. 271 § 49. Апостериорные оценки точности .............. 277 § 50. Некоторые замечания .................... 283 Приложение I. О распределении ошибок округления . . . 286 Приложение II. Решение больших задач линейной алгебры 293 Литература ............................. 301 Предметный указатель ....................... 302 »Там, где большинству алгориф- мистов-любителей кажется, что алгорифм готов для публикации, профессионал понимает, что тя- желая и утомительная работа только начинается» Дяс. Форсайт ПРЕДИСЛОВИЕ Сложная ли наука линейная алгебра? Начало пятидесятых годов. Первые выпуски студентов- вычислителей, перед которыми открывается увлекательный мир решения новых еще никем неизведанных проблем. И вдруг вместо заманчивой перспективы неожиданное пред- ложение—заняться созданием программного обеспечения электронных вычислительных машин (ЭВМ) для решения задач линейной алгебры. Особого восторга оно не вызвало. Легко понять, почему это произошло. Мы были вос- питаны в духе классических курсов, читаемых на мате- матическом факультете. Линейная алгебра была препод- несена нам столь четко и ясно, что не оставалось ни- каких сомнений в том, что все основные задачи, рас- сматриваемые етой областью математики, полностью решены. В самом деле, теория определителей исчерпывающе отвечала на вопрос о том, когда существует решение си- стемы линейных алгебраических уравнений, а правило Крамера указывало его явный вид. Все спектральные задачи сводились в основном к двум задачам — опре- делению корней алгебраического многочлена и решению систем уравнений. Более того, в нашем, арсенале были та- кие «эффективные» численные-методы, как метод Гаусса, метод Данилевского и др. Эти методы вроде бы позволяли решать соответствующие задачи линейной алгебры во всей их полноте. Поэтому порученная нам работа первона- 6 ПРЕДИСЛОВИЕ чально воспринималась как чисто механический процесс по переводу огромного количества известных к тому времени вычислительных алгорифмов с общепринятого языка ма- тематических формул на язык команд ЭВМ. Действительность оказалась значительно сложнее, Лишь после многих неудач и ошибок мы стали понимать, что рядом с классической линейной алгеброй не только существует, но и успешно развивается совсем «другая» линейная алгебра, о которой почти ничего не говорилось ни в основных, ни даже в специальных курсах. Эта линей- ная алгебра была тесно связана со многими областями математики, уходила своими корнями в самые разнооб- разные приложения, заставляла учитывать особенности ЭВМ и языков программирования, требовала решения новых системных задач и никак не согласовывалась с ши- роко распространенным мнением о всемогуществе ЭВМ. Называлась она. «вычислительной», хотя данный термин далеко не полностью отражал содержание этой «другой» линейной алгебры и нередко низводил ее до уровня жон- глирования математическими преобразованиями. Сравнительная простота теории линейной алгебры и кажущаяся эффективность существовавших численных методов долгое время держали нас в своем плену. К со- жалению, многие математики и сейчас попадают под их успокоительное обаяние, не замечая всей сложности, ко- торая характерна для задач алгебры. Причина подобного положения кроется, на наш взгляд, в основах обучения этой науке, в методике -препода- вания, в содержании обязательных и специальных курсов линейной алгебры, читаемых в вузах. Вычислительная алгебра сделала за последние пятнадцать лет громадный скачок вперед и является одним из самых развитых на- правлений численного анализа. Содержание же лекций, как правило, слабо отражает достигнутый прогресс и по-прежнему ведется в духе изложения различных фак- ПРЕДИСЛОВИЕ 7 тов типа теорем существования без учета проблем вычис- лений. Знакомство с линейной алгеброй в вузе начинается для студентов-вычислителей с первых же лекций. Поэто- му от того, что и как читается в теоретической части этого курса, во многом зависит формирование основы будущего восприятия всей вычислительной математики. Нельзя не признавать красоту и изящество теории, построенной г;а понятиях линейной зависимости, базиса, определителя и т. п. Но все практические вычисления, связанные с ними, весьма неустойчивы. Поэтому методы исследования, используемые в теоретической части курса, оказываются не очень полезными при непосредственном их применении для конструирования численных методов и часто приводят просто к неправильному пониманию вычислительной стороны дела. Однако наличие подобных фактов в действительности предоставляет лектору исключительно благоприятную возможность формирования научных взглядов тех сту- дентов, для которых вычислительная математика как предмет должна занять существенное место в образова- нии. Конечно, реализация этой возможности требует из- менения всего теоретического курса линейной алгебры, но выигрыш от такого изменения может быть очень боль- шим. Численные методы алгебры станут естественной частью общего курса, не нужно будет дополнительно тратить драгоценные часы занятий на изложение их основ и, что самое главное, можно будет легко показать студенту громадную практическую значимость курса ли- нейной алгебры во всей ее полноте. Отсутствие органического единства в методике изло- жения теоретической и практической частей курса линей- ной алгебры не позволяет достичь должного эффекта в обучении студентов-вычислителей. В. полной мере мы ощутили это на себе в первые годы работы. С тех пор 8 ПРЕДИСЛОВИЕ прошло много лет. Но описанная выше ситуация с уди- вительным постоянством повторяется снова и снова. К со- жалению, и сейчас молодой специалист в области вычис- лительной математики нередко оказывается совершенно беспомощным, встретившись с необходимостью грамотно решить систему линейных алгебраических уравнений, не говоря уже о проблеме собственных значений. Все эти причины побудили нас предпринять попытку подготовить ряд взаимосвязанных и построенных на еди- ной основе учебных пособий по линейной алгебре, содер- жащих необходимый минимум теоретических знаний, без которого невозможно воспринять огромное вычислитель- ное богатство, и отражающих современные проблемы в области конструирования численных методов. Настоящее учебное пособие представляет собой третью книгу из этой серии после теоретического курса [1] и задачника [4]. Посвящено оно вычислительным основам линейной алгебры. Основной замысел книги был навеян теми трудностя- ми, с которыми приходится сталкиваться, изучая числен- ные методы линейной алгебры. Уже при первом знаком- стве с существующей литературой, например, по библио- графическому указателю [8], возникает недоуменный вопрос: «Почему решению в общем-то небольшого числа различных задач линейной алгебры посвящено такое ог- . ромное количество работ?» Нельзя на него дать однознач- ный ответ, ибо это явление обусловлено многими причи- нами. Для линейной алгебры характерна исключительная широта ее приложений. Учет конкретных особенностей задачи приводит к появлению новых модификаций чис- ленных методов, а желание решить данную задачу как можно лучше значительно увеличивает их число. И вот здесь, на наш взгляд, кроется одна из основных причин обилия публикаций. ПРЕДИСЛОВИЙ 9 Что значит решить задачу лучше? Если задача ре- шается на ЭВМ, то типичной ситуацией для линейной алгебры является использование стандартных программ. По крайней мере, так должно быть. Но для пользователя ЭВМ безразлично, какой из численных методов заложен в основу той или иной стандартной программы. Его интересуют, как правило, лишь три ее характеристики: время счета, объем.требуемой памяти ЭВМ и точность. Относительно легко сравнить методы по первым двум характеристикам. Что же касается точности, то это уже сделать значительно сложнее. Трудно даже ответить на вопрос, как сравнивать методы по точности. Именно этим обстоятельством и объясняется наличие большого числа работ, в которых либо' ничего не говорится о точности численных методов, либо доказательство преимущества одних из них сводится к неубедительным эмпирическим доводам и эмоциональным рассуждениям. Мы уже отмечали обманчивость простоты формулиро- вок задач линейной алгебры. Однако во всей полноте это постигается лишь тогда, когда проводится анализ влия- ния ошибок округления и возмущения входных данных на точность решения. Существенный прогресс в исследовании устойчивости численных методов произошел сравнительно недавно и связан с возникновением так называемого обратного ана- лиза ошибок. Основная идея этого анализа заключается в том, что реально вычисленное решение рассматривается как точное для той же задачи, но с возмущенными вход- ными данными. При этом само возмущение выбирается так, чтобы его действие оказалось эквивалентным сово- купному влиянию всех ошибок округления. Обратный анализ предложил идею, но не инструмент для изучения ошибок округления. Даже для самых про- стых алгорифмов исследование эквивалентных возмуще- ний остается тяжелой и утомительной работой, сопро- 10 ПРЕДИСЛОВИИ вождающейся выполнением большого числа весьма тонких выкладок. Тем не менее обратный анализ позволил оце- нить совместное влияние ошибок округления и ошибок входных данных на точность результатов и провести на этой основе сравнение численных методов между собой. Исследование лучших численных методов линейной алгебры показало удивительную малость их эквивалент- ных возмущений. Большая часть этих методов, предназ- наченных для решения систем уравнений и проблемы собственных значений, представлена в настоящей книге. В ней же описаны и многие вспомогательные алгебраи- ческие алгорифмы, позволяющие конструировать новые численные методы. Возможно, что после знакомства с этой книгой у чи- тателя появится желание использовать для решения своей задачи какой-либо иной метод. Прежде чем реализовать такое желание, стоит, выполнить для нового метода столь же тщательный анализ ошибок, какой проделан здесь для всех численных методов. Что же касается ответа на вопрос, поставленный в на- чале предисловия, то линейная алгебра действительно простая наука, если оставаться в пределах классических формулировок обязательного курса и не замечать тех сложных проблем, которые в ней существуют. А каково ваше мнение? В. В. Воеводин ГЛАВА I МАТЕМАТИЧЕСКИЕ ОСОБЕННОСТИ МАШИННОЙ АРИФМЕТИКИ • Современная вычислительная техника стала необхо- димым звеном выполнения самых различных научных исследований. Она позволяет автоматизировать сложней- шие вычислительные процессы и получать достаточно быстро и в нужной форме решение многих задач. Однако за всем этим в действительности скрыто преобразование огромной информации. Это преобразование может быть весьма сложным или совсем простым, но в конечном счете оно всегда сводится к выполнению последователь- ности простейших операций, описанных системой команд электронной вычислительной машины. Общение с вычислительной техникой на уровне системы команд не эффективно для подавляющего большинства пользователей ЭВМ. Поэтому оно осуществляется на уровне каких-либо специальных машинно-независимых языков типа алгола, фортрана и других. Такие языки содержа г многие математические символы, с помощью которых принято описывать арифметические операции над числовыми, данными. Однако это не означает, что арифметические операции на ЭВМ обладают теми же свойствами, что и математические операции. Машинная арифметика имеет свои характерные особен- ности. Правильно учитывая их, можно достичь высокой эффективности в решении задач на ЭВМ. Невнимание же к этим особенностям нередко приводит к ошибочным результатам. § 1, Позиционные системы счисления Общий эффект от решения задачи и даже возможность ее решения во многом определяется тем, как в действи- тельности выполняются операции над числами. А это )2 ОСОБЕННОСТИ МАШИННОЙ АРИФМЕТИКИ (ГЛ I в свою очередь зависит от принятой системы записи чисел или, как говорят, системы счисления. Наиболее совершенным принципом записи чисел являет- ся тот, на котором основана наша десятичная система счисления. Известно, что любое неотрицательное чис- ло х может быть представлено в виде степенного ряда х==ац- Ю^+йя-г Ю^+.-.+йо+а-г lO-'+a-a- Ю-^... где коэффициенты а, могут принимать значения 0, 1,2,..., .... 9. Перечислив подряд все -коэффициенты, указав положение запятой и приписав числу некоторый знак, мы приходим к следующей системе записи: х = ± a„a;„-i... OQ, a-iQ-a... Несмотря на кажущуюся простоту, такая система явилась продуктом длительного исторического развития. Известный французский математик и физик Лаплас писал: «Мысль выражать все числа девятью знаками, придавая им, кроме значения по форме, еще значение по месту, на- столько проста, что именно из-за этой простоты трудно понять, насколько она удивительна. Как нелегко прийти к этой методике, мы видим на примере величайших гениев греческой учености Архимеда и Аполлония; от которых эта мысль осталась скрытой». Создание современной цифровой вычислительной тех- ники не связано с какими-либо принципиально другими системами счисления. Запись чисел, с которыми опери- рует ЭВМ, основана на той же идее, что и десятичная система. В математическом плане основные изменения невелики и заключаются в следующем. Зафиксируем некоторое целое положительное число р > 1 и целые числа Ко, «i, .... Kp-i. Пусть любое не- отрицательное число х может быть представлено в виде ряда x=bnpn+b^p't-l+...+b,+b.,pl+b.,p-2+..., (1.1) где каждый из коэффициентов &, может принимать одно из значений Сто, «i,..., ctp-i. Снова перечислив подряд все коэффициенты, указав положение запятой и приписав числу некоторый знак, мы получим аналогичную запись: x=±b„b^...b»,b ib г... (1.2) § 1] ПОЗИЦИОННЫЕ СИСТЕМЫ -СЧИСЛЕНИЯ 13 Описанные системы счисления называются позицион- ными. Их название связано с тем, что роль, которую играет каждое число в записи (1.2), зависит от занимаемой им позиции. Отсчет позиции определяется положением за- пятой или, что то же самое, положением коэффициента &о. В литературе, связанной с вычислительной математикой, слово «позиция» чаще всего заменяется словом «разряд». Нумерация разрядов устанавливается в убывающем порядке подряд слева направо, причем первый разряд слева от запя- той имеет нулевой номер. Различаются разряды числа до запятой ц разряды после запятой. Число р называется основанием системы счисления, числа «о, a.i, ..., йр-i — базисными. Если используется система счисления с осно- ванием р, то правую-часть (1.2) называют р-ичной дробью. Дробь называется бесконечной, если в ее записи (1.2) имеется бесконечно много ненулевых коэффициентов, и ко- нечной в противном случае. Обычно в записи дроби (1.2) опускаются все первые и последние нулевые коэффициенты. Опускается и запятая, если все коэффициенты после нее являются-нулевыми. Выбор базисных чисел йо, osi, ..'., Kp_i определяется в основном требованиями удобства работы с веществен- ными числами в данной системе счисления. Не видно каких-либо особых преимуществ, которые дало бы исполь- зование базисных чисел, превосходящих по модулю осно- вание системы счисления. Поэтому мы будем считать, что \^и\<Р (1.3) для всех k. В современной цифровой вычислительной технике чаще всего используются системы счисления с базисными числами к/,==А. Теорема 1.1. Если. базисные числа образуют сово- купость О, 1, ..., р—\,пго любое вещественное число мо- жет быть представлено в виде р-ичной дроби (1.2). Доказательство. Покажем, что любое число х может быть представлено в виде ряда (1.1). Очевид- но, что достаточно рассмотреть лишь положительные числа х. Существует целое число п^ такое, что выполняются соотношения p'^^-six0. Находим далее целое число «з такое, что pn^f^x—ЬntPt^^ «2. Если х — Ь„,р"' — Ьп.р^ == О, то получение ряда (1.1) закончено. Поэтому снова рас- сматриваем случай х - Ьп,р^ - Ьп,р^ > 0. Продолжая этот процесс, получаем последовательность целых чисел rii> п^> Пз>... и чисел Ьп„ Ьп„ Ьп„ ..., выбираемых из совокупности 1, 2, ..., р—\. При этом, либо при некотором k х-^Ьп/^О, (1.5) i=i либо для всех k k p"^i-^р5, (2.2) либо Xs, либо Xs-^-p3, если \x—Xs\==n Р5- Замена числа х числом ^ является операцией округ- ления, причем лучшей во многих отношениях. Однако по сравнению с описанной ранее операцией она имеет два существенных недостатка. Как следует из (2.2), для ее реализации необходимо выполнять операцию сложения примерно в половине случаев. Так как округление числа осуществляется после каждой арифметической опе- рации, его реализация согласно (2.2) фактически приводит к замедлению работы арифметических устройств ЭВМ. Кроме этого, есть некоторая неоднозначность в выполнении данного варианта операции округления, связанная с третьим случаем из (2.2). Мы увидим в дальнейшем, что она не так уж безобидна. Естественным является желание объединить достоин- ^ ства обоих способов округления. Покажем, что этого 18 ОСОБЕННОСТИ МАШИННОЙ АРИФМЕТИКИ [ГЛ Г можно добиться путем использования специальных систем счисления. До сих пор мы предполагали, что в качестве базисных чисел р-ичной системы счисления используются числа О, 1, ..., р — 1. При этом оказалось, что лучший по точности способ округления чисел в таких системах не является самым простым в реализации и приводит к замедлению выполнения арифметических операций. Рассмотрим теперь р-ичные позиционные системы счисления с другими набо- рами базисных чисел. Прежде чем исследовать возмож- ность представления чисел в таких системах, посмотрим, чего можно достичь выбором базисных чисел. Пусть числа задаются в р-ичной системе с базисными числами Ok. Наилучшее приближение х;? для х и в этой системе дает оценку ошибки округления на классе веще- ственных чисел не лучше, чем (2.1). Поэтому, если мы найдем базисные числа такие, чтобы для всех чисел х выполнялось неравенство х,-х ^-^-Р5, (2.3) то эта система счисления во многих отношениях будет наилучшей. Именно, операция округления, заключающаяся в отбрасывании всех разрядов, начиная с s—1-го, будет не только самой простой, но и будет иметь наименьшую ошибку округления. Для выполнения условия (2.3) необходимо и достаточно, чтобы все числа вида О ...О й,_Д-2 ... были по модулю не более (1/2)?". Для этого в свою оче- редь необходимо и достаточно выполнение условия max \a, -^?-1. (2.4) 0<*<Р—1 ~ Любая р-ичная система счисления должна иметь р раз- личных базисных чисел. Но неравенство (2.4) имеет р различных целочисленных решений относительно cc,i, лишь при нечетном р, при этом сами кд определяются одно- значно. Именно, a,=(l+2fe-p)/2. Таким образом, если искомая система счисления суще- ствует, то она должна иметь нечетное основание и базисные § 2] ОКРУГЛЕНИЕ ЧИСЕЛ 19 числа—(1/2) (р— 1), ...,+(1/2) Co—I). Заметим, что в таких системах не нужно дополнительного изображения для знака числа, так как он учитывается в его цифровой части. Подобные системы счисления называются сокращенными. Теорема 2.1. Если р—нечетное и базисные числа образуют совокупность — (1 /2) (р — 1), ..., + (1 /2) (р — 1), то любое вещественное число х может быть представлено в виде ряда (1.1). Доказательство. Будем снова считать для опре- деленности, что число х положительное. Доказательство этой теоремы также основано на последовательном опре- делении всех коэффициентов ряда (1.1) и в идейном плане почти не отличается от доказательства теоремы 1.1. Един- ственное отличие состоит в следующем. При получении ряда (1.1) в теореме 1.1 конечные дроби всегда приближали число л; с .недостатком. Теперь же коэффициенты ряда (1.1) находятся из условия, чтобы соответствующая конечная дробь наилучшим образом приближала число л", т. е. допускается как приближение, с недостатком, так и с из- бытком. Мы не будем останавливаться на доказательстве бодее подробно, оставляя проведение его читателю в каче- стве упражнения. Среди сокращенных позиционных систем счисления простейшей является троичная сокращенная система. Как уже отмечалось, в современной вычислительной технике наиболее широко используется двоичная система. С точки зрения округления чисел этот выбор не является лучшим, так как, например: Троичная сокращенная система счисления обладает более простым способом наилучшего округления чисел и не при- водит к замедлению выполнения арифметических операций. В дальнейшем мы покажем и другие преимущества сокращенных систем счисления с точки зрения влияния ошибок округления на вычислительный процесс. УПРАЖНЕНИЯ Всюду рассматривается сокращенная система счисления. 1. Написать р-ичную дробь числа р. 2. Доказать, что знак числа совпадает со знаком первого разряда. 3. Доказать, что х^~>х, если первый из ненулевых отброшенных разрядов отрицательный, и х^<_х, если он положительный. 4. Доказать, что число (1/2) р5 не может быть задано конечной дробью.