Оценки погрешности методов Ланцоша и Арнольди в точной и машинной арифметике тема автореферата и диссертации по математике, 01.01.07 ВАК РФ

На правах рукописи

УМА-'

Книжнерман Леонид Аронович

Оценки погрешности методов Ланцоша и Арнольди в точной и машинной арифметике

01.01.07 — вычислительная математика

АВТОРЕФЕРАТ диссертации на соискание учёной степени доктора физико-математических наук

- 1 НОЯ 2012

Москва — 2012

005054213

Работа выполнена в Федеральном государственном бюджетном учреждении науки Институте вычислительной математики Российской академии наук

Официальные оппоненты:

доктор физико-математических наук, профессор Валерий Павлович Ильин (Федеральное государственное бюджетное учреждение науки Институт вычислительной математики и математической геофизики Сибирского отделения Российской академии наук, лаборатория вычислительной физики, г. н. с.)

доктор физико-математических наук Игорь Евгеньевич Капорин (Федеральное государственное бюджетное учреждение науки Вычислительный центр им. А. А. Дородницына Российской академии наук, с. н. с.)

доктор физико-математических наук Сергей Павлович Суетин (Федеральное государственное бюджетное учреждение науки Математический институт им. В. А. Стеклова Российской академии наук, отдел комплексного анализа, в. н. с.)

Ведущая организация: Федеральное государственное бюджетное учреждение науки Институт проблем управления им. В. А. Трапезникова Российской академии наук

Защита состоится 2012 г. в -^Ч часов на заседании дис-

сертационного совета Д 002.04^.01 в Федеральном государственном бюджетном учреждении науки Институте вычислительной математики Российской академии наук по адресу: 119333, г. Москва, ул. Губкина, д. 8.

С диссертацией можно ознакомиться в библиотеке Федерального государственного бюджетного учреждения науки Института вычислительной математики Российской академии наук.

Автореферат разослан

тт^о- о12 г-

Учёный секретарь

диссертационного совета Д 002.045.01 ^ у у

доктор физико-математических наук ^-/ы^О^ Г. А. Бочаров

Общая характеристика работы

Исторический обзор

Численное решение дифференциальных и интегральных уравнений после дискретизации часто сводится к выбору приближённого решения из подпространства Крылова1

К.т{А, <р) = зрап{Л V Л V, • • •, Ат~1<р}, (1)

где. Л — матрица размера п х п и р — вектор (далее предполагаемый нормированным). При непосредственной работе со степенным базисом

(2)

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

1. Осуществляется в неявном виде ортогонализация по Граму-Шми-дту базиса (2). Соответствующие методы в случаях эрмитовой и неэрмитовой матрицы А называются методами Ланцоша и Арнолъди. Проведение классической (эрмитовой) ортогонализации позволяет при исследовании этих методов использовать технику ортогональных многочленов и многочленов Фабера. Если А — не матрица, а ограниченный оператор в гильбертовом пространстве, то для получения оценок можно работать со спектрами и использовать теорию логарифмического потенциала. Вклад в обоснование именно этих методов призвана внести представленная диссертация.

2. Все другие случаи: строится квазиортогональный базис или пара базисов. Методы этого подсемейства семейства крыловских методов (неэрмитовы обобщения метода Ланцоша и др.) также используются на практике, но в них случаются неприятности, когда в знаменателях расчётных формул появляются маленькие числа или ноль. Обычно эти методы можно использовать, когда имеется хороший предобуславливатель, но и это не даёт гарантии сходимости. Путей получения общих априорных оценок для методов этого подсемейства не видно.

'А. Н. Крылов, Изв. АН СССР, VII сер., Отд. матем. и ест. наук, Кг 4 (1931), 491-539.

аТ. е. не требующих знания результатов вычислений, в отличие от апостериорных-, см. А. А. Амосов, Ю. А. Дубинский, Н. В. Копченова, Вычислительные методы для инженеров, М., Высшая школа, 1994: с. 61-62. Под погрешностью метода мы понимаем отклонение искомой величины от выдаваемого методом приближения (а не оценки промежуточных величин). См. с. 23-24 там же.

Метод Арнольди3 вычисления собственных пар неэрмитовых матриц и несамосопряжённых операторов стал востребованным в 1980-х годах благодаря4 работам Ю. Саада. Метод Арнольди применяется во многих задачах: собственно для вычисления спектра;5 8 как средство решения систем линейных алгебраических уравнений;7 8 для локализации спектра в итерационных методах решения линейных систем;9 10 при решении жёстких систем обыкновенных дифференциальных уравнений;11 12 для вычисления матричной экспоненты.13 14 15 16

Пусть А — ограниченный оператор в гильбертовом пространстве <р е Н, ||у>|| = 1. Процесс Арнольди в И с А и <р основан на ор-тогонализации по Граму-Шмидту степенного базиса подпространства Крылова (1). Первые т шагов процесса Арнольди можно выразить соотношением

AQ = QH + hm+hmqm+ie^, (3)

где Q = (qi,..., qrn) — набор первых т ортонормированных векторов Арнольди, верхняя хессенбергова т х m-матрица Н = (h¡j) содержит коэффициенты рекурсии, a есть j-й единичный вектор размерности т).

Числа Ритца (собственные значения Н), которые считают приближениями к собственным значениям А, не обязаны лежать в спектре S оператора А\ они находятся в «поле значений» (множестве значений отношения Рэлея) {{А"ф,ф) I гр 6 ||?/>|| = 1}. Замыкание поля значений будем обозначать К.

Предпринятые до наших работ попытки оценить скорость сходимости метода Арнольди для матриц не привели, на наш взгляд, к получению законченных естественных результатов. Расстояние от собственного вектора до крыловского подпространства (1) оценено Ю. Саадом17 через коэффициенты разложения <р по собственным векторам А и значе-

3W. Е. Arnoldi, Quart. Appl. Math., 9 (1951), 17-29.

4X. Д. Икрамов, Несимметричная проблема собственных значений, М., Наука, 1991: § 1С.

6А. Ruhe, Lect. Notes in Math., 973 (1982), 104-120.

eY. Saad, Linear Algebm and its Applic., 34 (1980), 269-295.

7Y. Saad, Matii. Сотр., 37 (1981), 105-126.

8Y. Saad, M. H. Schultz, SIAM J. Sei. Stat. Сотр., 7 (1986), 856-869.

9H. C. Elmau, Y. Saad, P. E. Saylor, SIAM J. Sei. Stat. Сотр., 7 (1986), 840-855.

10D. Ho, F. Chatelm, M. Beunanii, Math. Mod. and Numer. Anal, 24 (1990), 53-65.

llP. N. Brown, Y. Saad, SIAM J. Sei. Stat. Сотр., 11 (1990), 450-481.

12C. W. Gear, Y. Saad, SIAM J. Sei. Stat. Сотр., 4 (1983), 583-001.

>3E. Gallopulos, Y. Saad, SIAM J. Sei. Stat. Сотр., 13 (1992), 1236-1264.

14M. Hochbruck, C. Lubich, J. Numer. Anal., 34 (1997), 1911-1925.

l5M. Hochbruck, C. Lubich, H. Selhofer, SIAM J. Sei. Сотр., 19 (1998), 1552-1574.

16Y. Saad, SIAM J. Numer. Anal., 29 (1992), 209-228.

1TY. Saad, .Numerical methods for large eigenvalue problems, Manchester, Manchester Univ. Press, 1992: глава VI, § 7.

ния многочленов степени не выше то — 1 в точках спектра:

dist (zj,K.m(A,ip)) (4)

< Г 4' min max |р(А)|,

где orjt — коэффициенты разложения tp = Ylk=i ce^Zk входного вектора tp по набору нормированных собственных векторов z^ матрицы А. Коэффициенты oik могут быть велики при больших перекосах где-то в спектре А. Таким образом, плохая обусловленность -«неинтересной» моды может послужить причиной «развала» оценки для искомой собственной пары. Кроме того, оценить расстояние от собственного вектора до крыловского подпространства недостаточно: из малости этого расстояния не следует малость расстояния от искомого собственного вектора до какого-либо вектора Ритца. Стьюарт18 распространяет результат Саада на недиа-гонализируемые матрицы, но оценивает то же расстояние, и его оценка также может содержать необоснованно большие величины. Стьюартом и Джа19 погрешность аппроксимации собственного значения Xj на шаге т оценена по порядку корнем m-ой степени из расстояния от соответствующего собственного вектора z} до подпространства (1):

min — <4 (Ä "" (2|И| + -РШ """

(5)

где 5 = sin Z(zj,fCm(A, tp)). Глава 9 диссертации показывает, что в одной из типичных ситуаций указанное расстояние убывает со скоростью геометрической прогрессии при росте тп, а тогда из оценки (5) вряд ли можно извлечь что-либо полезное.

' Для сравнения скажем, что одна из наших простейших оценок погрешности решения спектральной задачи имеет следующий вид: если есть конечное количество s собственных значений Aj,..., As, не принадлежащих полю значений К сужения оператора А на подпространство, дополнительное к интересующим нас s модам (это условие даёт возможность вывести из оценок расстояний от нужных собственных векторов до крыловского подпространства оценки ошибок соответствующих векторов Ритца), то оценка погрешности вычисления Aj (1 < j < s) равна 0{р™тс), где числа 0 < pj < 1 определяются взаимным расположением A j и К, с — константа, зависящая от формы границы К, и под знаком О не скрыты никакие перекосы (перекосы влияют на К). В общем

18G. W. Stewart, Matrix algorithms. Volume II: eigensystems, Philadelphia, SIAM, 2001: гл. 4, теорема 3.10.

19Zh. Jia, G. W. Stewart, Math. Сотр., 70 (2001), 637-647: следствие 4.2.

случае оценка не может быть существенно улучшена. Из этого следует, что оценки (4) и (5), игнорирующие поле значений сужения оператора, в принципе не могут даже гарантировать факт сходимости при т —»■ оо для оператора или семейства матриц, удовлетворяющих условиям нашей оценки.

В тот же период Ю. Саадом20 была дана оценка погрешности решения системы линейных уравнений Аи = <р с помощью обобщённого метода минимальных невязок GMRES, использующего базис, построенный процессом Арнольда.21 В этой оценке также (как и в (4)) фигурируют коэффициенты разложения правой части <р по системе собственных векторов А.

В [4], [5, § 4], [18, § 5] автором получены оценки погрешности метода спектрального разложения Арнольди (МСРА), т. е. варианта метода Арнольди, предназначенного для вычисления произведения матричной (операторной) функции на вектор. Если функция / аналитична в окрестности S и окрестности Sp (Я) — спектра Я, то МСРА в качестве приближения к вектору

и = f(A)ip (6)

выдаёт вектор

ит = Qf(H)ei. (7)

Из рекурсии (3) следует, что МСРА точен для многочленов степени < т — 1, то есть

f(A)<p = Qf(H)elt fe C[i], deg/<m-l. (8)

Если / аналитична на всём компакте К, то её можно разложить там в ряд Фабера и благодаря (8) свести оценку погрешности МСРА к оценке остатка ряда Фабера. Использование поля значений А в статье [4] благодаря оценке нормы резольвенты вне К освободило нас от необходимости явным образом отражать в оценке погрешности метода перекосы в спектре. Случай функций, которые могут иметь особенности между несколькими «внешними» собственными значениями А, рассмотрен нами в работах [5] и [18]; там с аналогичной целью мы использовали поле значений сужения А на подходящее инвариантное подпространство. А в работе [19] мы получили результаты, учитывающие спектр оператора (здесь — не матрицы) Л; в частности, для FOM (метода полной ортогонализации вычисления A~lip) показана сходимость на некоторой последовательности шагов, если 0 лежит в неограниченной связной компоненте C\Sp(A), о чём не было речи в старой матричной теории.

2°Y. Saad, Iterative methods for sparse linear systems, Philadelphia, SIAM, 2003: предложение 6.32.

"Г. И. Марчук, Ю. А. Кузнецов, Докл. АН СССР, 181 (1968), № 6, 1331-1334.

В [5, § 2], [18, § 3] мы получили оценки для точности вычисления собственных пар. Первое, что хочется сделать для оценки качества выделения собственной пары (A, z), ■— воспользоваться формулой (8) и взять многочлен / степени не выше тп — 1, равный 1 в точке Л и принимающий как можно меньшие значения па остальной части спектра или на поле значений ограничения А на подходящее инвариантное подпространство («дополнительное» к Cz). Так и делается, но ситуация осложняется ненормальностью оператора А и ненормальностью Я. Часть оценок представляет собой обычные неравенства, но некоторые оценки получаются коши-адамаровского типа: оценивается верхний или нижний предел то-го корня ошибки в терминах функции Грина (с особенностью в бесконечности) поля значений или спектра. Нижним пределом приходится-довольствоваться, когда поле значений ограничения А на инвариантное подпространство содержит искомое собственное значение. Наши оценки для спектральной задачи также сформулированы в терминах упомянутых функций Грина для подходящего сужения оператора; соответственно, в доказательствах используются элементы теории потенциала. То, что было сказано выше об отсутствии необходимости явно учитывать перекосы в спектре и об оценках, верных для подпоследовательностей, справедливо и здесь. Использование в формулировках функции Грина позволяет рассматривать произвольные, а не частные (отрезки, эллипсы) поля значений.

В статье [21] мы для частного случая f(A) = (zI—A)~l и д = 1 устгъ новили в терминах операторов жёсткие ограничения на Л, при выполнении которых есть толк от квадратуры Гаусса-Арнольди,22 23 под которой понимается приближённое равенство {f(A)tp,g(A)ip} « {f(H)ei,g(H)ei), где функции /ид аналитичны на спектрах Л и Я; до нас оставался открытым вопрос о том, быстрее ли сходится квадратура по сравнению с приближённым вычислением векторов f(Ä)<p и g(A)tp с помощью метода спектрального разложения Арнольда.

Метод Ланцоша24 теоретически представляет собой метод Арнольда, применённый к самосопряжённому оператору: А* = Л. В этом случае матрица Н коэффициентов рекурсии является вещественной симметричной трёхдиагональной, а многочлены Арнольди Q; (Qi(A)ip = qi+1) ортогональны на вещественной оси и удовлетворяют трёхчленному рекуррентному соотношению. Компакт К превращается в отрезок вещественной оси; ряд Фабера функции / на К становится сдвинутым рядом Фурье

22D. Calvetti, S.-M. Kim, L. Reichel, SIAM J. Matrix Anal. Appl., 26 (2005), 765-781.

23R. W. Freund, M. Hochbruck, In: Linear Algebra for Large Scale and Real-Time Applications¡ M. S. Moonen et al. (eds), Kluwer Academic Press, Amsterdam, 1993, 377-380.

24C. Lanczos, J. Res. Nat. Bur. Standards: Sect. В, 45 (1950), 225-280.

по многочленам Чебышёва. Эрмитовость матрицы H сильно облегчает построение теории в точной арифметике. Общая оценка погрешности метода спектрального разложения Ланцоша (MCPJI; эрмитов аналог MCP А) была получена (для точной арифметики) в [2]; она позволила оценить ошибку вычисления произвольного собственного значения методом Ланцоша. До наших работ были известны лишь теоремы Каниэля25 и Саада,26 уточнённые Пэйджем,27 которые касались аппроксимации со скоростью геометрической прогрессии нескольких собственных значений на одном краю спектра. Результаты Каниэля и Саада оценивают качество приближения j-го точного собственного значения А¡ j-ым же числом Ритца Ôj. Оценки эффективны, если число этих собственных значений существенно меньше, чем размерность подпространства Крылова т.

Используя равенство (8), В. JT. Друскин и автор получили априорные оценки ошибки вычисления не обязательно близкого к краю спектра собственного значения методом Ланцоша.

Исследование метода Ланцоша существенно усложняется при переходе к машинной арифметике (естественно, в случае матриц).28 29 30 Было выяснено, что векторы Ланцоша теряют ортогональность и даже становятся почти линейно зависимыми после того, как появится хорошее приближение к какому-либо собственному значению. Процесс Ланцоша с одной и той же парой (Л, <р) может по-разному идти на компьютерах разных типов или даже на одном компьютере при разных способах компиляции программы. Причина в том, что в процессе Ланцоша при стандартной реализации вследствие теоретической трёхдиагональности H новый вектор Ланцоша ортогонализуется только к двум предыдущим, а не ко всем, в отличие от метода Арнольди. Поведению процесса Ланцоша в условиях машинной арифметики посвящены основополагающие работы К. Пэйджа,31 32 33 в которых виртуозно вскрыт механизм появления неустойчивости в процессе Ланцоша и влияние неустойчивости на базовые матричные соотношения процесса. Оценки погрешности34 можно извлечь

"S. Kaniel, Math. Сотр., 20 (1966), 369-378.

2eY. Saad, SI AM J. Numer. Anal, 17 (1980), 687-70G.

"Б. Парлетт, Симметричная проблема собственных значений, М., Мир, 1983: § 12.4, (12.4.17).

2аС. К. Годунов, Г. П. Прокопов, Ж. вычиел. матем. и матем. физ., 10 (1970), 1180-1190.

29Б. Парлетт, Op. cit.: гл. 13.

30G. H. Golub, D. P. O'Leary, SIAM Review, 31 (1989), 50-102.

31C. C. Paige, The computation of eigenvalues and eigenvectors of very large sparse matrices, Ph. D. thesis, London, Univ. of London, 1971.

32C. C. Paige, J. Inst. Math. Applic., 18 (1976), 341-349.

33C. C. Paige, Linear Algebra and its Applic., 34 (1980), 235-258.

34Говорж о погрешности хрьшовских методов в машинной арифметике, мы имеем в виду суммарный эффект от остановки на конкретном шаге m и от ошибок округления.

из эстетичной работы Гринбаум,35 которая, используя работы Пэйджа, результат реального процесса Ланцоша на шаге тп представила как результат идеального процесса Ланиоша с матрицей размерности п + тп, но в упомянутой работе наложены очень жёсткие ограничения на параметры задачи и нельзя в общем случае получить оценку погрешности лучше корня четвёртой степени е1//4 из элементарной ошибки округления

Оценка погрешности МСРЛ в машинной арифметике с правильным порядком по е была дана в [3]. Для анализа влияния ошибок округления на процесс Ланцоша было построено чебышёвское рекуррентное соотношение с правой частью порядка ошибки округления; эта правая часть отвечает за неточность процесса. Получение оценки ошибки МСРЛ, в правой части которой слагаемое, отвечающее за ошибки округления, линейно по е (что выведено из результатов Пэйджа), свелось к доказательству устойчивости чебышёвской рекурсии.

В [3], [11] был также в принципе объяснён феномен Ланцоша — свойство процесса Ланцоша рано или поздно вырабатывать приближения к собственным значениям с почти машинной точностью. Попытка объяснить феномен Ланцоша была предпринята Каллэм и Уиллафби,38 однако в их работе используется ряд недоказанных предположений.37

Априорных оценок погрешности МСРЛ и МСРА до наших работ не было, если не считать оценки для решения систем линейных уравнений (f(A) = А"1).

Итак, до нас:

• спектральные «оценки» в неэрмитовом случае или содержали в левой части промежуточные величины (невязки, расстояния от собственного вектора до крыловского подпространства), а не истинные отклонения чисел и векторов Ритца от собственных значений и векторов соответственно, или в правой части необоснованно содержали потенциально большие величины (результаты Саада, Стьюарта и др.). Непосредственно применимые (а не содержащие оценки лишь промежуточных величин) априорные результаты о сходимости чисел Ритца к собственным значениям относились лишь к краю спектрального интервала в эрмитовом случае при точной арифметике (оценки Каниэля и Саада);

• из матричных функций исследовалась только функция ¡{А) = Л-1, соответствующая решению системы линейных уравнений. Не было средств теоретического рассмотрения даже такой важной и популярной функции, как матричная экспонента f(A) = exp(-tA), t > О, А > 0;

35А. Greenbaum, Linear Algebra and its Applic., 113 (1989), 7-63.

3eJ. Cullum, R. A. Willoughby, Linear Algebra and its Applic., 29 (1980), 63-90.

"X. Д. Икрамов, Итоги науки и техн.: Матем. анализ, 20 (1982), 179-260: § 9.

• априорные оценки для простого процесса Ланцоша и квадратуры Гаусса-Ланцоша в машинной арифметике имели в лучшем случае правую часть порядка п3тех/А (результат Гринбаум и следствия из него);

• Были лишь частные попытки объяснить явление адаптации методов к спектру. Оценки для квадратуры Гаусса-Арнольди ничего не говорили о её эффективности: не было ясно, имеет ли (и если имеет, то при каких условиях на спектр) применение квадратуры преимущество перед вычислением соответствующих матричных функций.

Цель работы:

построение общей теории сходимости методов Ланцоша и Арнольди.

Актуальность темы

Методы Ланцоша и Арнольди — классические методы вычислительной математики. Они широко используются,38 им уделяется много внимания. Однако имевшиеся до нас априорные результаты об их сходимости носили частный или незаконченный характер. Это и обусловило актуальность темы диссертации.

Научная новизна

В диссертации систематически изложена общая теория сходимости методов Ланцоша и Арнольди (в том числе в форме MCP Л и MCP А), ликвидирующая имевшиеся качественные пробелы и содержащая прямые априорные оценки погрешности.

Теоретическая и практическая ценность

Результаты диссертации дают пользователям (математикам и программистам):

• принципиальную уверенность в том, что обсуждаемые методы при установленных условиях работают. В частности, прояснены классы матричных функций, к вычислению которых можно применять методы спектрального разложения Ланцоша и Арнольди;

• возможность априори оценивать объём вычислительной работы, достаточный для решения конкретной задачи, а также, зная вид аналитической зависимости границы ошибки от номера шага т, апостериорно прогнозировать фактическую ошибку путём подбора параметров (curve fitting);

38На момент написания данного данного предложения поисковая система Google выдаёт 79 тысяч ссылок на запрос "Lanczos method" и 16 тысяч — на запрос "Arnoldi method".

• инструмент для исследования методов в условиях машинной арифметики: возмущённые чебышёвские и фаберовские рекуррентные соотношения.

Ссылки на описания некоторых приложений, к которым имеет отношение автор, даны в диссертации. Приложений, возможно, и не было бы, если б не была развита теория.

Методы исследования

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

Положения, выносимые на защиту:

• Доказано, что при вычислении методом Арнольди «крайних» собственных пар сходимость имеет место со скоростью геометрической прогрессии, показатель которой выражается через значение функции Грина дополнения к полю значений сужения оператора на инвариантное подпространство, порождённое «неинтересной» частью спектра. Предложен метод спектрального разложения Арнольди для вычисления умноженной на вектор операторной функции; дана оценка погрешности этого метода в терминах коэффициентов ряда Фабера вычисляемой функции, построенного на поле значений.

• Исследован метод спектрального разложения Ланцоша для вычисления умноженной на вектор функции от эрмитовой матрицы. Для точной-арифметики дана оценка погрешности метода в терминах коэффициентов ряда вычисляемой функции по многочленам Чебышёва, смещённым на спектральный интервал. Оценена погрешность вычисления произвольно расположенного в спектральном интервале собственного значения методом Ланцоша; в частности, показано, что числа Ритца сходятся к заданному собственному значению со скоростью геометрической прогрессии, показатель которой определяется отделённостью собственного значения в спектре.

39Автор полагает, что он был среди тех, кто в начале 1990-х годов стал активно внедрять использование поля значений, функции Грина и многочленов Фабера в исследование сходимости крыловских методов. Из «соратников» отметим М. Айерманна, О. Неванлинну и Л. Н. Трефевена с соавторами (список не претендует на полноту). Сейчас этими понятиями удивишь не очень многих вычислительных линейных алгебраистов: они (понятия) начали входить в учебники, а тогда их использование встречало заметное сопротивление.

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

• В терминах линейных ограниченных операторов показано, почему на практике обсуждаемые методы работают лучше, .чем гарантировано предыдущими оценками (адаптация к спектру). Адаптивные оценки погрешности сформулированы в терминах; функции Грина неограниченной компоненты дополнения к операторному спектру. Получены результаты о сходимости на подпоследовательностях шагов процесса Ланцоша или Арнольди в ситуациях, когда процесс сходиться в обычном смысле не обязан. Проанализирована эффективность квадратуры Гаусса-Арнольди.

Обоснованность и достоверность результатов

Представленные в диссертации утверждения снабжены математическими доказательствами. Результаты численных экспериментов иллюстрируют многие из доказанных утверждений. Кроме того, справедливость положений представленной теории подтверждается многолетней практикой моделирования геофизических и иных процессов с помощью методов Ланцоша и Арнольди.

Апробация работы

Результаты диссертации излагались на юбилейной ланцошевской конференции (Рэйли, США, 1993), конференции по вычислительной линейной алгебре (Миловы, Чехия, 1997), конференции 81АМ-САММ-2006 (Дюссельдорф, Германия), на кафедре вычислительной математики Технического горного университета Фрайберга (Германия), в ВЦ РАН им. А. А. Дородницына, в ИПУ РАН им. В. А. Трапезникова, по нескольку раз обсуждались на семинарах в Институте вычислительной математики РАН и Институте прикладной математики РАН им. М. В. Келдыша (в том числе на семинаре им. К. И. Бабенко).

Публикации

Основные результаты диссертации изложены в статьях [2, 3, 4, 5, 6, 10, 11, 13, 17, 18, 21], опубликованных в рецензируемых научных журналах из списка, рекомендованного ВАК. Другие результаты диссертации опубликованы в работах [7] и |19] (последняя реферирована в гег^гаШа«). Всего по теме диссертации, включая аннотации конференционных докладов и препринты, автором опубликована 21 работа, из них 12 — в рецензируемых научных изданиях.

В статьях, написанных совместно с В. Л. Друскиным, лично диссертант: использовал разложение в ряды Фурье-Чебышёва и специальные функции (В. Л. Друскин оценивал коэффициенты нужных рядов Фурье-Чебышёва с помощью сведения к оценке погрешности разностных схем' для обыкновенных дифференциальных уравнений40); провёл выкладки, связанные с многочленами Чебышёва; разработал и использовал технику возмущённых чебышёвских рекуррентных соотношений. Включённые в диссертацию результаты, полученные совместно с В. Л. Друскиным и Э. Гринбаум, принадлежат их авторам в равной степени.

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

Ссылки на статьи автора по теме диссертации можно найти в нескольких книгах41 42 43 44 45 46 47 и в статьях разных авторов.

Благодарности

Автор благодарит: В. Л. Друскина — за многолетнюю совместную деятельность по исследованию и практическому применению метода Лан-цоша, а также за полезные обсуждения личных статей автора; Э. Гринбаум — за совместное написание нескольких статей по теории метода Ланцоша, а также за полезные обсуждения личных статей автора; за полезные обсуждения диссертации в целом или части отражённых в ней статей и за библиографическую поддержку — М. Айерманна, С. К. Асва-

дурова, Б. Бекерманна, А. Б. Богатырёва, Дж. Голуба , С. А. Горейнова,

40V. Druskin, L. Knizhnerman, Radio Science, 29 (1994), 937-953: § 5.

4,C. К. Годунов, Лекции no современным аспектам линейной алгебры, Новосибирск, Научная книга (ИДМИ), 2002: предисловие,

"G. H. Golub, G. Meurant, Matrices, momenU and quadrature with applications, Princeton and Oxford, Princeton Univ. Press, 2009: § 7.5, 8.0.

"A. Greenbaum, Iterative methods for solving linear systems, Philadelphia, SIAM, 1997: гл. 4.

"N. J. Higham, Functions of matrices. Theory and computation, Philadelphia, SIAM, 2008: § 13.5-13.6.

45G. Meurant, The Lanczos and Congugate Gradient algorithms: From theory to finite precision computations, Philadelphia, SIAM, 2006: § 3.11.

4eB. N. Parlett, The symmetric eigenvalue problem, Philadelphia, SIAM, 1998: § 13.11.

"L. N. Trefethen, M. Embree, Spectra and pseudospectra: the behavior of nonnormal matrices and operators, Princeton, Princeton Univ. Press, 2005: гл. VI, § 28.

С. Н. Давыдычеву, Н. JI. Замарашкина, X. Д. Икрамова, В. П. Ильина, Дж. Каллэм, А. В. Князева,, К. Любиха, Ю. М. Нечепуренко,

B. С. Рябенького, А. Скороходова, А. В. Собяшша, В. 3. Соколинского,

C. П. Суетина, К. Торреса-Вердина, Е. Е. Тыртьннникова, М. Хохбрук, О. Эрнста, участников перечисленных в пункте «Апробация» конференций и семинаров; Дэвида Бейли — за публикацию фортранного пакета повышенной точности48 в Интернете; Центральную геофизическую экспедицию и Schlumberger-Doll Research — за организационную помощь.

Структура и объём работы

Диссертация состоит из двенадцати глав, десять из которых содержат математические результаты и распределены между четырьмя частями. Объём диссертации — 255 с. Список литературы содержит 187 наименований.

Краткое содержание работы

Глава 1 «Введение» описывает задачи, даёт исторический обзор и содержит информацию, обязательную для диссертаций.

Часть 1 «Методы Ланцоша и спектрального разложения

ланцоша в точной арифметике: неадаптивные оценки» Содержит главы 2-3. Термин «неадаптивные» означает «не адаптированные к спектру», то есть зависящие лишь от спектрального интервала оператора в случае метода спектрального разложения Ланцоша и не зависящие от лакун в спектре, не связанных с отделённостью нужного собственного значения, в случае метода Ланцоша. Адаптивные оценки рассматриваются позже.

Глава 2 «Метод спектрального разложения Ланцоша в точной арифметике» содержит начальную информацию про МСРЛ. Пусть А — симметричная вещественная матрица размера пхп, / — функция, определённая на спектральном интервале А, <р € R". Обозначим через А<, г = 1,2,..., n, Aj+i > А», собственные значения А (А„ > Ai), а через — соответствующую ортонормированную систему собственных векторов. Пусть разложение <р по собственным векторам имеет вид

п ¿=1

Напомним содержание т шагов процесса Ланцоша приближённого вычисления спектра оператора А. В подпространстве Крылова (1)

48D. N. Bailey, ACM Trans. Math. Software, 21 (1995), 379-387.

строится базис 91,..., дт, получаемый в результате ортогонализации по Граму-Шмидту последовательности уэ, А<р,..., Ортогонализацию

можно провести с помощью трёхчленной рекурсии

= + г = 1,2,... ,тп — 1,

где /?о?о = 0, дх = <р и /?( > 0. Обозначим через Н трёхдиагональную

\

матрицу

Н =

(а 1 А

А ot2 ßi

\ ßm-l Olm)

через в{, Si = (sij,..., smi)T, i = 1,2, ...,m, — собственные значения и соответствующие нормированные собственные векторы Н. Положим Q = {qi...qm), yi = Qsi. Тогда () — приближённые собственные пары оператора А, получаемые методом Ритца на Кт(А, ф). Положим

А„ + Ai -о = --—in

-А, д(х) = /

+ Ai — (А„ — Ai)a;

An — Ai А n — Ai

где In — единичная матрица порядка п. Рассмотрим ряд Фурье-Чебышёва функции д:

д{ х) = Y,9kTk{x).

(10)

k=0

Мы оцениваем скорость сходимости МСРЛ через скорость сходимости ряда (10).

Теорема 1. Пусть ряд (10) абсолютно сходится на отрезке [—1,1]. Тогда для вектора (7) справедлива оценка погрешности

II" - ""mit < 2 \дк\ < +оо.

(11)

к=т

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

Предложение 1. Пусть Аа > 0, f(A) = ехр(-М), Ь > 0, а = При т < а справедлива оценка погрешности

U-«m <2

Vi+ °(7)

у/Б

т

■ехр

-f+0

2а

т

Предложение 2. Пусть Ai > 0, f(A) = cos , г = x

f = ^jp — 1. Если 0 < £ < 1, то имеет место оценка погрешности

1 + 0(1 + ^)ехр{-х [igf! + 0(a]}

у/тт</Ат2 -х2 уД

Введём в рассмотрение обратную функцию Жуковского Ф(гу) — w + Vw2 -1.

Предложение 3. Пусть Ai > 0, f(A) = ехр (-zVA), z > 0. Справедлива оценка погрешности

П-1 А„ + А;

||и-ит||<4Ф(с)-т[1-Ф(с)-1,] \ с =

А„ — Ai

Предложение 4. Пусть А, > 0, Ап < 1, /(А) = А3, б € М\{0}. Справедлива оценка погрешности

11«-«»||< ^[и-оН+оГ-^-^ехрГ-^ + О^)

тп V 7Г \s-mj в

Глава 3 «Метод Ланцоша в точной арифметике» содержит доказательство аппроксимационных оценок для собственного значения, не зависящих от его номера в спектре. Первая из оценок не учитывает отделённость желаемого собственного значения в спектре, а вторая учитывает.

Предложение 5. Предположим, что ||Л|| < 1, Аг = 0, <рг ф 0 и т > 3. Тогда

GM™'1

l<i<nJ^ ~ 2-* ^

(\ l°g(A~2)

ф — 2J, то правая часть (12) близка к —^— Определим величину

6 = min ( max |/(А)| f(X) € ВД, deg> < т - 1, /(Аг)

^ASIAbAr-ilUlAr+i^n]

где 1<г<п, ^г>0ит>3. Теорема 2. Если ||Л|| < 1 и

6< *

l+VT^' 16

то выполняется оценка

шш |Лг-6>,1< f^H

, . я^ -7-(13)

1 <j<m J ' ' л v '

Правая часть оценки (13) убывает со скоростью геометрической прогрессии благодаря оптимизации многочлена /,49 Б0, например, по Золотарёву.

Часть 2 «Методы Ланцоша и спектрального разложения Ланцоша в машинной арифметике» состоит из глав 4-6 и также содержит неадаптивные оценки. Оценки из глав 4-5 являются аналогами оценок из глав 2-3 соответственно для машинной арифметики. Оценки. главы 6 для квадратуры Гаусса-Ланцоша являются аналогами хорошо известных свойств квадратур типа Гаусса.

Глава 4 «Метод спектрального разложения Ланцоша в машинной арифметике» напоминает начальную информацию о машинной арифметике и содержит машинно-арифметические оценки погрешности МСРЛ. Буквой е обозначается элементарная машинная ошибка округления. Будем придерживаться относительно машинной арифметики предположений из книги Парлетта.51 Обозначения Я и т. п. будут относиться к соответствующим реально вычисленным объектам. Введём также следующие обозначения: с\ — максимальное количество ненулевых элементов в строках А, е0 = 2(п + 4)е, ег = ^7 + е> е2 = у/2тах(б£о,£'1), V = т2,5[|Л||е2- Пэйджем показано, что при условиях

7л[6(п + 4)г + £]] <1, (п + 4)с<^ (14)

справедливо неравенство

\1-г]<61<\п + т), г = 1,2,..., тп. (15)

Как принято в теории метода Ланцоша, мы считаем, что все действия с трёхдиагональной матрицей Н, включая решение спектральной задачи, производятся точно.

Сделаем спектральное линейное преобразование, переводящее отрезок [А1 — т], А1 + т)}, заведомо содержащий, в силу (15), собственные значения Н, в отрезок [—1,1]. А именно, положим

'С*п + ^1) - (Ап - Ах + 2г})х'

9(х) = f

49В. И. Лебедев, Ж. вычисл. матем. и матем. физ., 9 (1969), 1247-1252.

60Е. М. Никишин, В. Н. Сорокин, Рациональные аппроксимации и ортогональность, М., Наука, 1988: гл. 2, доказательство утверждения 8.3].

51Б. Парлетт, Op. cit.: § 2.4.

Введём в рассмотрение ряд Фурье-Чебышёва

оо

g{x) = Y,9kTk{x), -1<®<1. (16)

к=О

Теорема 3. Пусть функция f определена на интервале [Ai —7/, А„+ 77] и ряд Фурье-Чебышёва (16) абсолютно сходится на отрезке [—1,1]. Тогда вектор ит определён корректно и при выполнении (14) справедлива оценка погрешности

[2 imii ig °°

||W - «rnll < ^-М^-гли^)^^^ + -щу/™ Е Ы- (17)

к—т .

Член оценки (17), отвечающий за машинную арифметику, имеет порядок машинной ошибки округления.

Пример запуска MCPJI. В Горном техническом университете Фрайберга (Германия) одна из научных программ, написанная автором совместно с В. JI. Друскиным, вычисляющая матричную функцию f(A) = A~l ехр(—tA), t > 0, была запущена с одними и теми же входными данными на персональных компьютерах, работавших под управлением операционных систем Windows и Linux. Трансляция была осуществлена фортранными компиляторами фирмы Intel: автором для Windows и в ИВМ РАН — для Linux. Таблица 1 показывает, что два процесса Ланцо-ша в конце концов выдали приближённое решение, совпадающее в пяти десятичных знаках. При этом коэффициенты ланцошевской рекурсии и значения приближённого решения по ходу процесса (от момента потери ортогональности векторами Ланцоша и до достижения сходимости заданного уровня) различались весьма значительно.

Глава 5 «Феномен Ланцоша и расположение чисел Ритца» содержит машинно- арифметические оценки погрешности процесса Ланцоша. Рассмотрим процесс Ланцоша с начальным вектором (9) и будем интересоваться качеством приближений к фиксированному собственному значению Аг матрицы А, 1 < г < п. Без ограничения общности можно считать, что ||Л|| < 0.9 и Аг = 0.

Предложение 6. Пусть ||Л|| < 0.9, Аг = 0, <рг ф 0 в (9) и т > 3. Предположим, что выполнены условия Пэйджа (14), а также условия

„ 25 - , 2^1 /0.19. .

9т £2 < 1, n»aei<-W—Ы.

Тогда справедлива оценка

(18)

min = min — Ar| < -

lg«<m' 1 l<t<m 1 - 2

\y/6\Vr\)

Таблица 1: Результаты вычисления матричной функции на двух компьютерах, работающих под управлением операционных систем Windows и Linux, j — номер шага процессов Ланпоша, aj и /3,- — коэффициенты ланцошевских рекурсий, Uj — значения приближённого решения в отдельном (контрольном) узле пространственной разностной сетки.__

3 Linux Windows

ßj aJ Uj ßl aj Uj

50 -0.15702E-14 -0.15702E-14

100 998.23 20434. -0.16543E-09 998.23 20434. -0.16543E-09

150 -0.15208E-07 -0.15208E-07

200 1987.9 2817.1 -0.79739E-07 1985.8 2813.8 -0.79744E-07

250 -0.62672E-06 0.63888E-06

300 5083.4 2987.7 -0.31798E-05 9944.3 9050.3 -0.32142E-05

350 -0.77681E-05 -0.76203E-05

400 10448. 12052. -0.11058E-04 8183.5 5388.9 -0.11047E-04

450 -0.11686E-04 -0.11686E-04

500 11734. 11345. -0.11666E-04 9273.5 5455.4 -0.11666E-04

Если

log , « 1,

т — 1 т/б\(р., то правая часть оценки (18) примерно равна

1 , 10 Jm

— log ——. т ч/б|<Рг{

Таким образом, оценка имеет по т. порядок

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

Теорема 4. Пусть 1 < г < п, <рг > 0, j = тшх<*<п, ¡¿г |АГ - А,-| > 0, т > 3, /(X) G Ж[Х\, deg/ < т - 3, /(Ar) = 1,

npU Х € ^ ~Т}'Хп + 7]}< ~ " при X 6 [Al-7?,ЛП + Г7]\]ЛГ_1 + i,Ar+i

и (помимо (Ц)) выполнены условия Б < с^рг и тег < с37||Л[|-1, где И = у/тпц + Д, 5 = с4т3£Г2||Л||, Д = т3е2||Л||(А„ - А^-1, и где с2 и Сз — достаточно малые, а с4 — достаточно большая положительная

константа. Тогда справедлива оценка

ЦЛЦ-1 min \0i - Ar|

1<»<т

< шах Упф? (m^^lHI^-1 + D) , + ^ЦЛЦт"1] .

При достаточно большом допустимом числе шагов оценка теоремы 4 имеет порядок машинной ошибки округления.

Когда матрица А является невырожденной, но не положительно определённой, то матрица Н — Нт, порождённая т шагами процесса Лан-цоша, может быть вырожденной. В этом случае аппроксиманта MCPJI к решению системы линейных уравнений Аи = <р не определена на т-ом шаге процесса. Мы показываем, однако, что по меньшей мере одна из двух последовательных трёхдиагональных матриц, Нт или #т_ 1} имеет спектр, чьё относительное расстояние до нуля больше, чем примерно нормированный квадрат расстояния от Sp(yl) до нуля.

Предложение 7. Положим

d(z) = ,_nün \z - А,-1 = dist(z, Sp(vl)).

Пусть в — собственное значение Нт, ц — собственное значение Нт-\, ||Л|| = 1 и ни в, ни [i не слишком близко к спектру А, именно,

min{d(0),d(At)} > [("г + I)3 + л/Зт2] е2.

Тогда

max{|0|, |М|} >-Щ*-а*1 ,щ

~ d(0) + 6.25/Зт + V(6-25/?m)2 + 12.5/3^(0)+^' V '

где а < 12.5 -I- 0(e). При условии

90.625 [(2п + 6)е + т(3г0 + <ñ)] + 12.5\/m£i < 8.5

неравенство (19) может быть заменено на

тах{|0|, И} > ® - + ОИ.

Попутно в этой главе получено обобщение теоремы Стилтьеса о еле] баум.54

разделении,52 53 улучшающее аналогичный результат из работы Грин-

"Т. J. Stieltjes, Annales scientifiques de E.N.S. 3e série, 1 (1884), 409-426.

53E. M. Никишин, В. H. Сорокин, Op. cit.: гл. 2, § 8, упражнение 6.

54A. Greenbaum, Op. cit.

Глава 6 «Гауссова квадратурная формула, порождаемая простым процессом Ланцоша, и её приложениям посвящена квадратурной формуле типа Гаусса (Гаусса-Ланцоша, КФГЛ), порождаемой т первыми шагами простого (т. е. без реортогонализации) процесса Ланцоша. Под этим понимается квадратурная формула (для матриц это, точнее, формула приближённого суммирования)

т

(f(A)qi,qi) « (f(H)euei) =

i=1

где / — функция, определённая и гладкая на отрезке [Aj - 77, An + г]] и где первое скалярное произведение берётся в R", а второе — в Rm. Нарушение ортогональности векторов Ланцоша в простом процессе Ланцоша не позволяет просто обосновать законность использования КФГЛ, порождаемой m первыми шагами простого процесса Ланцоша, и побуждает использовать вычислительно дорогую реортогонализацию, что нежелательно. Это порождает задачу обоснования КФГЛ в условиях машинной арифметики.

Теорема 5. Пусть а = ||А||/(АП - Ax). Пусть функции fug определены на [Ai -rj, А„ + 77] и разлагаются в ряды

- 2А

k=Q ч-v »1 + 21?

сходящиеся абсолютно на указанном отрезке. Тогда справедливы неравенства

\(ПА)дид(А)Я1) - (/(Я)е1|5(Я)е1)| <0(т3ае2) £ \/ш\ + 2 £ \/ш\

к+1<2т-2 к+1>2т-1

и

^НАЫ'-ЫЮе^КО^ае^ £ |М|+2 £ |/,/,|.

к+1<2т-2 к+1>2т-1

Следствие 1. В условиях теоремы 5 верна оценка

2т—2 оо

1</(ЛЬ,?1) - {}{Н)еие1)I < О (т*ае2) £ |Д| +2 £ |Д|. (20)

0 к=2т—1

Бесконечный ряд в оценке (20) для КФГЛ начинается примерно с вдвое большего номера, чем в оценке (17) теоремы 3 для МСРЛ.

Рис. 1: Пример вычисления величины (А íq^.,q\)■ Кривая 1 — график 1о£10 кривая 2 — 1сщ10 91) - {Я-^.е!)!, кривая 3 — 1с^10 - (С}Я^е!,^)!.

В конце главы описываются приложения доказанных оценок: конструктивные — к методу спектрального разложения Ланцоша и решению некорректных задач с помощью вариационной регуляризации и теоретические — к феномену Ланцоша. В иллюстративных целях рассматривается пример со сравнением вычисления ф) с помощью КФГЛ и непосредственно с помощью МСРЛ. Матрица А положительно определена и имеет хорошо отделённое собственное значение. Приведены результат ты счёта для чётных т (для нечётных т картинка несколько сдвинута, но качественно не отличается). Рисунок 1 показывает, что величина (Я-1е 1,в1) стабильно сходится к решению, в то время как ((¿Н^ех, дх) делает «подскоки», связанные, очевидно, с моментами потери ортогональности, происходящей при сходимости очередного дубля к хорошо отделённому собственному значению.

Часть 3 «Методы Арнольди и спектрального разложения Арнольди: неадаптивные оценки» содержит главы 7-9. «Неадаптивность» здесь понимается в том смысле, что оценки сформулированы в терминах поля значений оператора или сужения оператора на подходящее инвариантное подпространство, а не в терминах спектра. Здесь мы получаем аналоги результатов части 1 для несамосопряжённого случая; вместо многочленов Чебышёва и рядов Фурье-Чебышёва появляются многочлены и ряды Фабера; появляются также функции Грина, в терминах которых можно оценить поведение нужных экстремальных многочленов.

В главе 7 «Метод спектрального разложения Арнольд и: общие оценки» оценивается погрешность вычисления вектора (6) с помощью метода Арнольди в несимметричном случае в условиях точной и машинной арифметики. Вместо разложения / в ряд Фурье-Чебышёва, применявшегося в симметричном случае, здесь используется разложение в ряд Фабера55 на поле значений.

Процесс Арнольди в компьютерной арифметике с реортогонализа-цией может быть выражен формулами AQ - QH = /im+i,mîm+iem + F и QTQ — I = G, где под Q, H и др. понимаются реально вычисленные на ЭВМ объекты и нормы п х n-матрицы F и m х т-матрицы G — умеренные кратные величины е. Мы даём оценки погрешности в терминах величин ||F||, ||G|| и ||QTgm+i||, не вдаваясь в стандартную для вычислительной линейной алгебры задачу оценивания упомянутых норм. (Величины ||F|| и max(||G||, ||QTÇm+i||) не превосходят произведений, соответственно, £||А|| и е на одночлены от m, п с небольшими показателями степеней и коэффициентами.) Положим ri = ||F|| + max(||G||, ||Qrgm+i||)IHI- Предположим, что выполняются условия ||G|| < 0.5, \\F\\ < ||Л||, \\QTqm+i\\ < 1, ||?m+i|| < 2, |/im+i,m| < 2||Л||. Они естественны и следуют из оценок точности вычисления скалярных произведений, линейных комбинаций векторов, умножения матрицы на вектор и т. п.56 57 при наложении на m, п и е ограничений типа полиномиальных неравенств.

Свяжем с матрицей А множество К её значений отношения Рэлея. Будем рассматривать задачу вычисления (6) в предположении аналитичности / на К. Обозначим через ip функцию, конформно и биективно отображающую внешность единичного круга {ш G С | |и>| > 1} на дополнение к К в расширенной комплексной плоскости С при условиях 1р(оо) = оо и V'(°°) > °> через Г — границу К, а через Тр, р > 1, — линию уровня {ip(w) 11TiiI = р]. Через к е N, обозначаем многочлены Фабера для К — целые части k-oïi степени Пусть

ос

/(А) = Х>Ф*(А) (21)

fc^o

— ряд Фабера функции / на К.

Теорема 6. Если функция f аполитична в канонической области Gr, ограниченной кривой Гд с параметром R > 1+7?1/", то справедлива

55П. К. Суетин, Ряды по многочленам Фабера, М., Наука, 1984.

5вБ. Парлетт, Op. cit.: гл. 2.

"G. Н. Wilkinson, Rounding errors in algebraic processes, N. Y., Prentico-Hall, 1964.

оценка погрешности

,2(1+2.5

max

0<1<т-1

max (1 +-+ ч»аЛ (22)

;/<т-1 \ тах(/, 1) J

00 / i \k

где константа а > 1 и константа под знаком < зависят от К.

Следствие 2. Если е = 0 (точная арифметика) и функция / аполитична на К, то справедлива оценка погрешности

где константа а и константа под знаком < зависят от К.

С помощью недавнего результата Бекерманна58 можно убрать степенной множитель из правой части (23).

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

В методе Арнольди ортогональность векторов Арнольди устойчиво поддерживается, . поэтому вопросы машинной арифметики гораздо менее актуальны для метода Арнольди, чем для метода Ланцоша, частным случаем которого метод Арнольди является лишь в точной арифметике. Эта глава — последняя, где фигурирует машинная арифметика.

Глава 8 «Оценки погрешности метода Арнольди и метода спектрального разложения Арнольди: случай нормальной матрицы и крайних изолированных собственных значений» начинает изучение влияния дискретности края спектра на поведение процесса Арнольди. В предыдущей главе метод Арнольди, применённый к матрице А € Мп(Ш) и вектору <р е R", исследован как средство приближённого вычисления векторов вида (6), где / — функция, аналитическая на множестве отношений Рэлея матрицы А. В настоящей главе изучается сходимость в методе Арнольди «крайних» собственных значений, не принадлежащих выпуклой оболочке К остальных собственных значений, и

»8B. Beckermann, С. R. Acad. Sei. Paris: Ser. I, 340 (2005), 855-860.

59Хотя устойчивость вычисления конкретного базиса не гарантируется. См., например, A. Malyshev, М. Sadkane, Linear Algebra and its Applic., 371 (2003), 315-331; J.-F. Carpraux, S. K. Codunov, S. V. Kuznetsov, Linear Algebra and its Applic., 248 (1996), 137-160.

oo

(23)

k=m

соответствующих собственных векторов. Даны аналоги теорем Каниэ-ля и Саада из теории метода Ланцоша, сформулированные в терминах функции Грина (с особенностью в бесконечности) дополнения к К в расширенной комплексной плоскости. Затем исследуется метод Арнольди как метод вычисления (6) при ослабленном по сравнению с предыдущей главой ограничении на функцию /: ей разрешается иметь особенности между крайними собственными значениями А (/ должна быть аналитич-на в области, содержащей К и искомые собственные значения). Показано, что при больших т оценённая скорость сходимости МСРА не зависит от выделенных собственных значений.

Глава 9 «Оценки погрешности метода Арнольди и метода спектрального разложения Арнольди: случай крайних изолированных собственных значений» содержит результаты, аналогичные результатам предыдущей главы, но не использующие предположения о нормальности. В терминах функции Грина дополнения к полю значений К сужения оператора на подпространство, дополнительное к искомым собственным векторам, даны более слабые (что обусловлено отсутствием теоремы о невязке с учётом отделённости) аналоги теорем Каниэля и Саада и последнего результата предыдущей главы.

Пусть U — гильбертово пространство, А: Н->Н — ограниченный линейный оператор, ip — элемент П с ||у>|| = 1. Мы рассмотрим процесс Арнольди с An ¡р. Пусть АЬ...,Л, (s > 1) — различные собственные значения А, £ = [Щ^Й " V)] <Р € Н, Q — замыкание линейного подпространства span {£, А£, ...}. Будучи замыканием инвариантного подпространства, Q само инвариантно, ц мы можем определить ограничение В = А\д: Q G и замыкание К поля значений В.

Мы исследуем сходимость собственных значений Аь ..., Xs в предположениях, что Ai <£ К,..., As $ К и что ^ ф О,..., ipa ф 0. Пусть Ф — функция, конформно и биективно отображающая дополнение к единичному кругу {w € С | Н > 1} на дополнение к К и С при условиях Ф(оо) = оо, Ф'(оо) > 0, Ф — функция, обратная к Ф, U = Ф(А,). Определим величины =| U |-щ TnCs/ifi, где константа с5 > 1 зависит от формы К.

Теорема 7. Метод Арнольди, применённый к оператору А и вектору ц>, зат > s шагов выработает такие приближённые собственные пары (01, yi),..., (вт, ут), что при подходящей нумерации и нормализации верны оценки

К-0г<5(, (г = 1,...,s).

Обозначим через Фк (к 6 N) многочлены Фабера, построенные на

К, через Гр — изолинии функции Грина Гр = { Ф(ги) | w € С, |u>| = р}, р > 1. По функции /, аналитической на К, можно построить ряд Фабера / = YX= о fk^k, сходящийся на К.

Теорема 8. Пусть функция f аполитична в области, содержащей внутренность контура Гд (R > I) и окрестности точек Ai,...,Ae. Тогда метод Арнольди, применённый к А и <р, за достаточно большое число шагов т выработает корректно определённую аппроксиманту ит к и, такую что

Пт ||u-um||1/m < Я"1.

т-+оо

Часть 4 «Методы Ланцоша, Арнольди, спектральцого разложения Ланцоша и спектрального разложения Арнольди: адаптивные оценки» содержит главы 10-11. Гл. 10 объясняет, почему результаты применения методов Ланцоша и Арнольди обычно лучше, чем обещано оценками из предыдущих глав. Вводится понятие адаптации обоих методов к спектру оператора (который в данном случае нельзя считать матрицей), то есть зависимости качества результатов от спектра. Приводятся количественные выражения утверждений о том, что лакуны в спектре ускоряют сходимость и что сходимость в значительной степени обусловлена спектром S, а не полем значений К. Гл. 11 в терминах операторов обсуждает качество квадратуры Гаусса-Арнольди для функции ((zl — Á)~vip, ¡p) и связь между методом Арнольди и рациональной аппроксимацией типа Паде функций марковского типа.

Глава 10 «Об адаптации методов Ланцоша и Арнольди к спектру, или почему два этих метода так мощны» трактует адаптацию изучаемых методов к операторному спектру. Простое вычисление начального куска ряда Фурье-Чебышёва (10) и ряда Фабера (21) с подставленной туда матрицей даёт несколько лучшие оценки, чем (11) и (23). Причина успешного использования МСРЛ/МСРА заключается в так называемой адаптации обоих методов к спектру, то есть в зависимости их эффективности от спектра А, теоретический учёт лакун в котором способен приводить к лучшим оценкам, чем (11) и (23). Данная глава посвящена как раз этому аспекту поведения МСРЛ/МСРА и собственно методов Ланцоша и Арнольди как методов вычисления спектра. По причине нерелевантности в ненормальном случае матричного спектра работа ведётся с операторами.

Пусть А — ограниченный оператор в гильбертовом пространстве Tí, Ч> &Н, S = Sp(A), D — связная компонента дополнения С\5, содержащая бесконечность, Е = С\D Э S, capí Е ■— логарифмическая ёмкость. Если caplS > 0, то существует обобщённая функция Грина g для Z);

в противном случае считаем, что g = +00. Обозначим через невязку FOM на к-ом шаге процесса Арнольди с А и <р (если приближённое решение MCP А не существует, мы полагаем = +оо).

Предложение 8. Если 0 £ Е, то для FOM справедлива оценка

lim ||^|Г/т<е-*(°>. (24)

m—>00

Адаптивный характер (24), так же, как и оценок ошибки из последующих утверждений, очевиден благодаря монотонности функции Грина60 по области: чем меньше множество Е, тем больше (нестрого) показатель д{0) и, значит, тем лучше сходимость. Обсуждаемая монотонность является строгой для регулярных областей. Так как К выпукло и содержит S, справедливо включение ECK. Если дк обозначает функцию Грина для К и Е ^ К, то д(0) > дк(0), что и выражает преимущество (24) над неадаптивной оценкой (23) ввиду свойств сходимости рядов Фабера.

Предложение 9. Пусть А = А* и 0 ^ S. Тогда

ïïmmintlIr^lMI^IO^^e-^), (25)

где g — обобщённая функция Грина спектра S.

Заметим, что в отличие от (24) оценка (25) гарантирует хорошее поведение FOM для последовательности значений m, имеющей положительную нижнюю натуральную плотность (> 0.5).

Отметим, что на практике метод Ланцоша используется для решения знакоопределённых и знаконеопределённых ССЛАУ в течение длительного времени.61 62 Однако при обосновании использования метода в случае знаконеопределённых операторов были трудности, связанные с проблемой устойчивости.63 Утверждая, что по крайней мере один из каждых двух последовательных ланцошевских шагов «хороший», предложение 9 как раз обосновывает применение MCPJI в этом случае (несмотря на возможные «выскоки»).

Определим «критический уровень функции Грина» Rq = inf| R >

1

К С Gr I, где G я — каноническая область для компакта Е.

60J. L. Walsh, Interpolation and approximation by rational functions in the complex domain, AMS Colloquium Publications, 20, Providence, RI, AMS, 1969: § 4.1.

6,B. N. Parlett, Linear Algebra and iU Applic., 29 (1980), 323-346.

62H. D. Simon, Math. Comp., 42 (1984), 115-142.

63C. C. Paige, B. N. Parlett, H. A. van der Vorst, Numer. Linear Algebra with Applic., 2 (1995), 115-133.

Теорема 9. Пусть R > Ro, а функция f аналитична в Gr и непрерывна на Gr. Тогда ошибка МСРА удовлетворяет неравенству

Пт Ilu - ит[|1/т < Л-1. (26)

т—«о

Если Е ф К, то (26) экспоненциально сильнее, чем (11) и (23).

Пусть 0 — изолированное простое собственное значение A, z — соответствующий собственный вектор, £ = Aip, Q — замыкание линейного подпространства span{£,Л£,Л2£, ...}, В = A\g: Q —» Q. В следующей теореме спектр S, множество Е и функция д относятся к оператору В (не А).

Теорема 10. Метод Арнольди, применённый к onepamöpy А и вектору <р, на т-ом шаге выработает приближённую собственную пару (9т,Ут), такую что при подходящей нормализации ут верны оценки

Иш \вт\х/т < е-а(0), Um ||z - ут\\^т < (27)

т—юо т—>оо

Заметим, что добавление конечного числа точек к спектру не меняет обобщённых изолиний64 функции Грина. Поэтому здесь нет необходимости совместно рассматривать несколько собственных значений, в отличие от теоремы 7.

Одна и та же последовательность значений т (которая может быть редкой) даёт два нижних предела в заключениях (27) теоремы 10.

Теорема 11. Пусть А = А* и 0 — изолированный элемент S (и, следовательно, собственное значение оператора А). Метод Ланцоша, применённый к паре (Л, <р), за т > 2 шагов выработает такие пары Ритца (вт,ут) с подходящим образом нормализованными векторами Ритца ут, что справедливы оценки

Ш \9т\1/т < е-®<°\ ЕЕ [min (|0m_!|, |0m|)]1/m <

m—>00 m—юо

BE [min (||Z - ym^\\, \\z - 3/m||)]1/m <

m->oo

где z — нормированный собственный вектор А, соответствующий собственному значению 0, g — обобщённая функция Грина для S и д(0) следует понимать как существующий предел <?(0) = Нтл_>о, л/о д(А) > 0.

MJ. L. Walsh, Op. cit.

Если ограничиться главной, экспоненциальной частью оценок, то можно сказать, что теорема 11 обобщает теоремы Каниэля и Саада в следующих двух направлениях. Во-первых, интересные собственные значения не обязаны лежать на краю спектра или быть соседними. Во-вторых, дополнительные (т. е. не связанные с отделённостью желаемых собственных значений) зазоры в спектре ускоряют сходимость, что как раз и выражает адаптацию. Отметим, что рассуждение на тему адаптации оценок Каниэля и Саада к спектру имеется в уже упомянутой книге Парлетта (§ 12.4).

Глава 11 «Квадратура Гаусса-Арнольди для функции {(zl— -A)-1V> V5) и Падe-подобная рациональная аппроксимация функций марковского типа» посвящена квадратуре Гаусса-Арнольди и связанных вопросам. Пусть А — ограниченный оператор в гильбертовом пространстве % и ip — нормированный вектор из И. Функция fm(z) = ((zl - Н)~1еие1) (z Sp(Я)) является [т — 1/ш]-аппроксимантой типа Паде в бесконечности с полюсами в точках 0¿ (1 < i < т) к функции f(z) = ((zl - A)~lip,ip) (z g S = Sp(j4)). Аппроксимация f(z) и fm(z) обсуждается в нескольких работах.65 66 Она является частным случаем квадратуры Гаусса-Арнольди. Пусть П — неограниченная связная компонента C\S, Е = C\f2 Э S.

Если оператор А нормален, то он обладает спектральным разложением,67 а пара (Л, <р) — спектральной мерой ц с носителем 5; имеем

S

Поскольку в терминах обобщённой функции Грина д области П мы имеем

Jm diet [(zl - Л)" V, fCm(A, <p)\1/m < e~*x\ z e П,

то два следующих утверждения являются критериями нетривиальности квадратуры Гаусса-Арнольди (мы считаем квадратуру нетривиальной, если с ростом т её ошибка уменьшается существенно быстрее, чем величина II(zl - - Q(zl - H)_1ei||).

Предложение 10. Пусть оператор А нормален и спектральная мера ß пары (A,ip) регулярна,68 Если существует щакая точка z € Q,

65D. Calvetti, S.-M. Kim, L. Reichel, SI AM J. Matrix Anal. Appl., 26 (2005), 765-781.

66R. W. FYeund, M. Hochbruck, In: Linear Algebra for Large Scale and Real-Time Applications, M. S. Moonen et al. (eds), Kluwer Academic Press, Amsterdam, 1993, 377-380.

бГУ. Рудин, Функциональный анализ, M., Мир, 1975: гл. 12.

МН. Stahl, V. Totik, General Orthogonal Polynomials, Encyclopedia of Math. Appl., 43, Cambridge, Cambridge Univ. Press, 1992: гл. 3.

что Ишт-юо |/(z) —/m(z)|1/,m < е g(z\ то у множества Е и, следовательно, у спектра S нет внутренних точек.

Теорема 12. Если спектр S нормального оператора А является частью аналитической кривой без самопересечений и не разделяет комплексную плоскость, то справедливы оценки

Ito I f(z) - fm(z)\l'm < e-»W, z 6 C\ Co(S),

m-tec

Um I f(z) - fm(z)|1/m < z € П = С\E,

m—>00

где Co(5) — выпуклая оболочка S.

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

Для марковской функции ß(z) = / где ц — положительная мера с компактным носителем бесконечной мощности S = Supp/i, лежащим в R, сходимость диагональных аппроксимант Паде гарантируется теоремой Маркова и её уточнениями;69 70 71 см. также препринт автора72 по поводу сходимости в лакунах спектра S по меньшей мере на каждом втором шаге ассоциированного процесса Ланцоша.

В работе Гончара73 (см. также работы Рахманова74 75) исследована сходимость поддиагональных аппроксимант Паде к функции ß + г, где г (Е C(i), г(оо) =0, — «рациональное возмущение» Д, при условии /¿'(i) > 0 почти всюду на выпуклом замыкании S (из чего следует, что S — отрезок); даны оценки сходимости полюсов аппроксимант к полюсам г. В статье Рахманова76 рассмотрен случай носителя меры — объединения нескольких отрезков вещественной прямой и вещественного рационального возмущения г 6 R(i). В статье Гончара и С. П. Суетина77 получены аналогичные результаты для комплекснозначных мер некоторых типов на вещественной прямой.

63А. Markoff, Acta Math., 19 (1895), 93-104.

70Е. М. Никишин, В. Н. Сорокин, Op. cit.: гл. 2, § 6.

71Н. Stahl, V. Totik, Op. cit.: гл. 6.

™L. Knizhnerman, Schlumberger-Doll Research, Res. Note EMG-001-95-12 (1995), 18 p.

"А. А. Гончар, Матем. сб., 97 (139) (1975), 607-629.

74E. А. Рахманов, Матем. сб., 103 (145) (1977), 237-252.

75Е. А. Рахманов, Матем. сб., 118 (160) (1982), 104-117.

78Е. А. Рахманов, Матем. сб., 104 (146) (1977), 271-291.

"А. А. Гончар, С. П. Суетин, Современные проблемы математики, 5 (2004), 3-67.

В статье С. П. Суетина78 отмечено, что если 5 не лежит на прямой, то предельные точки полюсов диагональных аппроксимант Паде могут быть всюду плотны в С. Эта и родственные трудности в теории аппроксимации Паде возникают из-за квази-, а не эрмитовой ортогональности знаменателей аппроксимант. Тот факт, что соответствующие квазиортогональные многочлены «хорошо себя ведут», приходится выводить из алгебраической специфики случая.79 80

Мы представляем рациональное возмущение с простыми полюсами в виде, удобном для проведения процесса Арнольди. На рациональное возмущение налагается условие вещественности суммарного вычета и его положительности.

Предложение 11. Пусть э > 1, Х1.....А, — попарно различные

комплексные числа, сь... ,с, — ненулевые комплексные числа, такие что 0 < сх-Ь • -+С, 6 К. Существует такая трёхдиагоналъная матрица Т е С", что <(*/ - Т)-^,в1> = £и * 6 С\{Ах,..., А5}.

Оставшиеся утверждения этой главы показывают, с какой скоростью^ метод Арнольди локализует полюсы рационального возмущения

с> = 1, Л > 0, функции марковского типа (28). Мы приведём здесь одно из них.

Матрицу Аг из условия предложения 12 можно найти с помощью предложения 11.

Предложение 12. Пусть ц — положительная мера с компактным носителем Б С С, А2 — оператор умножения на независимую переменную в Ь2^(3), К = 00(5) — замыкание поля значений А2, в > 1, числа Ах,..., Ав е П попарно различны, си...,с3 — ненулевые комплексные числа, сх + • • • + са = 1, в. > 0, матрица Аг е <С5Х* такова, что

3

<(*/ - Лх^ех, ех) = £ -2-, * 6 С\{Ах,..., Аа},

~ нормированный собственный вектор Аь отвечающий собственному значению А,, Н = С'фадЯ), А = ЛхфЛ2 - ортогональная прямая сумма пространств и соответствующая прямая сумма операторов. Процесс Арнольди с оператором А в И и начальным вектором

С. П. Суетин, Успехи матем. паук, 57 (2002), 45-142: пункт 3 введения. 79С. П. Суетин, Матем. сб. 191 (2000), 81-114. 80С. П. Суетин, Матем. сб. 194 (2003), 63-92.

где 1 — постоянная функция «единица» из L,2tlA{S), зат> s шагов выработает такие пары Ритца (dj, yj), j — 1,... ,т, что при подходящей их нумерации верны следующие оценки. Если при данном j (1 < j < s) Аj ^ К, то

lim |А, - 0,|1/m < e-s(ai)-min,<t<.s{At) g^ | _ ,1/m <

m-юо J т—*оо J

В общем случае (Аj g Е) lim |Аj - Öj|1/m < e-9(Ai)-mm1<l!<s9(At)) _ < е-5(Л,)_ (29)

m—ЮО 771—Ю0

Здесь g — обобщённая функция Грина области Q с особенностью в бесконечности. -

В гл. 12 (заключительной) сформулированы краткие выводы, перечислены нерешённые задачи, упомянуты подходы других авторов и высказаны благодарности.

Машинной арифметике посвящены главы 4-6 и частично 7. Результаты глав 10 и 11 применимы только к ограниченным операторам в бесконечномерных гильбертовых пространствах. Результаты, касающиеся машинной арифметики, применимы только к матрицам. Остальные результаты применимы как к операторам, так и к матрицам или классам матриц неограниченной размерности, причём пространство обычно можно считать комплексным, даже если приведённая формулировка вещественна.

Публикации автора по теме диссертации

[1] Друскин В. П., Книжнерман JI. А. Использование операторных рядов по ортогональным многочленам при вычислении функций от самосопряжённых операторов и обоснование феномена Ланцоша. Деп. в ВИНИТИ. 1987. № 1535-В87. 47 с.

[2] Друскин В. Л., Книжнерман JI. А. Два полиномиальных метода вычисления функций от симметричных матриц // Ж. вычисл. матем. и матем. физ. 1989. Т. 29. № 12. С, 1763-1775.

[3] Друскин В. JI., Книжнерман JI. А. Оценки ошибок в простом процессе Ланцоша при вычислении функций от симметричных матриц и собственных значений // Ж. вычисл. матем. и матем. физ. 1991. Т. 31. № 7. С. 970-983.

[4] Книжнерман JI. А. Вычисление функций от несимметричных матриц с помощью метода Арнольди // Ж. вычисл. матем. и матем. физ. 1991. Т. 31. № 1. С. 5-16.

[5] Книжнерман JI. А. Оценки погрешности метода Арнольди: случай нормальной матрицы // Ж. вычисл. матем. и матем. физ. 1992. Т. 32. № 9. С. 1347-1360.

[6] Druskin V., Knizhnerman L. The Lanczos optimization of a splitting-up method to solve homogeneous evolutionary equations //J. Comput. Appl. Math. 1992. V. 42. № 2. P. 221-231.

[7] Druskin V., Knizhnerman L. Evaluation for Krylov subspace approximation to internal eigenvalues of large symmetric matrices and bounded self-adjoint operators with continuous spectrum // Schlumberger-Doll Research. 1992. Res. note EMG-001. 20 p.

[8] Druskin V., Knizhnerman L. Error bounds for Lanczos method to compute matrix-vector functions // Cornelius Lanczos International Centenary Conference. Conference Program. December 12-17, 1993. North Carolina State University, Raleigh, North Carolina. P. 85.

[9] Knizhnerman L. The quality of approximations to a separated eigenvalue and location of values" // Cornelius Lanczos International Centenary Conference. Conference Program. December 12-17, 1993. North Carolina State University, Raleigh, North Carolina. P. 90.

[10] Druskin V., Knizhnerman L. Krylov subspace approximation of eigenpairs and matrix functions in exact and computer arithmetic // Numer. Linear Algebra with Appl. 1995. V. 2. № 3. P. 205-217.

[11] Книжнерман Л. А. Качество аппроксимаций к хорошо отделённому собственному значению и расположение «чисел Ритца» в простом процессе Ланцоша // Ж. вычисл. матем. и матем. физ. 1995. Т. 35. № 10. С. 1459 -1475.

[12] Knizhnerman L. On adaptation of the Lanczos method to the spectrum // Schlumberger-Doll Research. 1995. Res. Note EMG-001-95-12. 18 p.

[13] Книжнерман Л. А. Простой процесс Ланцоша: оценки погрешности гауссовой квадратурной формулы и их приложения // Ж. вычисл. матем. и матем. физ. 1996. Т. 36. N- 11. С. 5-19.

[14] Knizhnerman L. On adaptation of the Arnoldi method to the spectrum 11 Schlumberger-Doll Research. 1996. Report EMG-001-96-03. 13 p.

[15] Knizhnerman L. On adaptation of the Lanczos and Arnoldi methods to the spectrum // Workshop on Iterative Methods and Parallel Computing, June 16-21, 1997, Milovy, Czech Republic. Abstracts.

[16] Druskin V., Greenbaum A., Knizhnerman L. Using nonorthogonal Lanczos vectors in the computation of matrix functions // SIAM J. Sci. Сотр. 1998. V. 19. № 1. P. 38-54.

[17] Greenbaum A., Druskin V. L., Knizhnerman L. A. On solving indefinite symmetric lineax systems by means of the Lanczos method // Ж. вычисл. матем. и матем. физ. 1999. Т. 39. № 3. С. 371-377.

[18] Knizhnerman L. Error bounds for the Arnoldi method: a set of extreme eigenpairs // Linear Algebra and its Applic. 1999. V. 296. P. 191-211.

[19] Knizhnerman L. Adaptation of the Lanczos and Arnoldi methods to the spectrum, or why the two Krylov subspace methods are powerful // Чебышёвский сборник. 2002. V. 3. № 2. P. 141-164.

[20] Knizhnerman L. The quality of the Gauss-Arnoldi quadrature for {(zl - A)~ V, ip) and application to Pade-like approximation, In: SIAM-GAMM ALA 2006, Dusseldorf. Abstracts of the Applied Linear Algebra Conference SIAM-GAMM, 2006.

[21] Книжнерман JI. А. Квадратура Гаусса-Арнольди для функции ((zI—A)~ltp, ip) и Паде-подобная рациональная аппроксимация функций марковского типа // Матем. сб. 2008. Т. 199. № 2. С. 27-48.

Отпечатано в ООО «Техноком» ИНН/КПП 7724609946/772401001 115477, г.Москва, ул.Бехтерева, д. 13, корп.З Принято к печати 7.08.2012г.

Отпечатано 7.08.2012г. Тираж 120 шт. Заказ 399/03

Краткое описание проблемы и основных результатов

1 Введение

1.1 Семейство крыловских методов

1.2 Метод Арнольди. Роль теоретического использования поля значений и операторного спектра.

1.3 Метод Ланцоша. Роль теоретического использования че-бышёвских рекуррентных соотношений.

1.4 Структура и обзор содержания диссертации.

Часть I Методы Ланцоша и спектрального разложения Ланцоша в точной арифметике: неадаптивные оценки

2 Метод спектрального разложения Ланцоша в точной арифметике

2.1 Введение

2.2 Метод операторных рядов Чебышёва.

2.3 Метод спектрального разложения Ланцоша: оценка погрешности в терминах коэффициентов смещённого ряда Чебышёва вычисляемой функции.

3.2 Оценка качества аппроксимации собственного значения, не учитывающая отделённость.61

3.3 Оценка качества аппроксимации собственного значения, учитывающая отделённость .62

3.4 Краткие выводы.64

Часть II Методы Ланцоша и спектрального разложения

Ланцоша в машинной арифметике 66

4 Метод спектрального разложения Ланцоша в машинной арифметике 66

4.1 Введение .66

4.2 Начальные сведения о простом процессе Ланцоша.67

4.3 Возмущённые чебышёвские рекуррентные соотношения . . 68

4.4 Оценка погрешности метода спектрального разложения Ланцоша в машинной арифметике.70

4.5 Пример запуска.73

4.6 Краткие выводы.74

5 Феномен Ланцоша и расположение чисел Ритца 76

5.1 Введение.76

5.2 Феномен Ланцоша без учёта отделённости.77

5.3 Кластеризация чисел Ритца.80

5.4 Феномен Ланцоша с учётом отделённости.83

5.5 О числе элементов промежуточного кластера.89

5.6 Числа Ритца на последовательных шагах Ланцоша.93

5.7 Результаты численных экспериментов.97

5.8 Конкретизированная форма теоремы 5.1.100

5.9 Краткие выводы.101

6 Гауссова квадратурная формула, порождаемая простым процессом Ланцоша, и её приложения 103

6.1 Введение .103

6.2 Скалярные произведения матричных многочленов Чебышёва105

6.3 Оценки погрешности гауссовой квадратурной формулы и родственные утверждения.109

6.4 Возможные приложения.113

6.4.1 Вычисление (f(A)(p,tp).113

6.4.2 Решение некорректных задач с помощью вариационной регуляризации .116

6.4.3 Феномен Ланцоша.117

6.5 Краткие выводы.118

Часть III Методы Арнольди и спектрального разложения

Арнольди: неадаптивные оценки 119

7 Метод спектрального разложения Арнольди: общие оценки 119

7.1 Введение .119

7.2 Описание метода.120

7.3 Разложение в ряд Фабера.122

7.4 Возмущённые фаберовские рекуррентные соотношения . . .124

7.5 Теорема об оценке погрешности.126

7.6 Пример: системы линейных уравнений.130

7.7 Пример: транспонированные жордановы клетки.132

7.8 Краткие выводы.134

8 Оценки погрешности метода Арнольди и метода спектрального разложения Арнольди: случай нормальной матрицы и крайних изолированных собственных значений 135

8.1 Введение .135

8.2 Аппроксимация крайних собственных пар нормальной матрицы .136

8.3 Невозможность существенного улучшения оценок предложения 8.1 в неэрмитовом случае.143

8.4 Вычисление матричной функции с особенностями, лежащими среди крайних собственных значений.146

8.5 Краткие выводы.153

9 Оценки погрешности метода Арнольди и метода спектрального разложения Арнольди: случай крайних изолированных собственных значений 155

9.1 Введение .155

9.2 Постановка задачи и вспомогательные результаты.156

9.3 Аппроксимация собственных пар с внешними собственными значениями.160

9.4 Резольвенты и сопутствующие матрицы.165

9.5 Вычисление матричной функции.170

9.6 Численные примеры.172

9.7 Краткие выводы.175

Часть IV Методы Ланцоша, Арнольди, спектрального разложения Ланцоша и спектрального разложения Арнольди: адаптивные оценки 178

10 Об адаптации методов Ланцоша и Арнольди к спектру, или почему два этих метода так мощны 178

10.1 Введение .178

10.2 Решение линейных уравнений — общий случай.180

10.3 Решение линейных уравнений — самосопряжённый случай 184

10.4 Оценка погрешности методов спектрального разложения Ланцоша и Арнольди.186

10.5 Вычисление собственной пары — общий случай .189

10.6 Вычисление собственной пары — самосопряжённый случай .192

10.7 Результаты численных экспериментов.196

10.7.1 Общий случай.196

10.7.2 Самосопряжённый случай .197

10.8 Краткие выводы.202

11 Квадратура Гаусса—Арнольди для функции ((zl—А)~гср, tp) и Паде-подобная рациональная аппроксимация функций марковского типа 203

11.1 Введение .203

11.1.1 Аппроксимация Паде функций марковского типа. Выделение полюсов.203

11.1.2 Связь квадратуры Гаусса-Арнольди и аппроксимации типа Паде с полюсами в числах Ритца .204

11.1.3 Структура оставшейся части главы .205

11.2 Квадратура Гаусса-Арнольди: случай нормального оператора .206

11.3 Квадратура Гаусса-Арнольди: случай ненормального оператора .211

11.4 Представление суммы конечного числа простых полюсов с вещественным положительным суммарным вычетом в виде, соответствующем квадратуре Гаусса-Арнольди.213

11.5 Комплексное рациональное возмущение функции марковского типа: выделение простых полюсов с помощью метода Арнольди.215

11.6 Численный пример с кратным полюсом.220

11.7 Краткие выводы.228

12 Заключение, задачи, другие подходы и благодарности 229

12.1 Заключение.229

12.2 Задачи.230

12.3 Другие подходы и родственные результаты других авторов 231

12.3.1 Подход Э. Гринбаум к теории простого процесса Ланцоша.231

12.3.2 Методы Ланцоша и Арнольди и теория ортогональных многочленов.232

12.3.3 Псевдоспектр (спектральный портрет матрицы) . . 233

12.3.4 Условные минимизационные задачи теории потенциала и адаптация к спектру.233

12.3.5 Методы с рестартами.234

12.3.6 Итерации подпространства.234

12.3.7 Крыловские процессы с преобразованным оператором235

12.4 Благодарности.235

Литература 237

Обозначение Смысл

Скалярное произведение

То же, что О

Обращение <

Равенство по определению; сравнение по модулю аддитивной подгруппы

Конец доказательства

Таблица 1: Обозначения (математические символы).

Сокращение Расшифровка

МЛИ Метод локальных итераций

МОРЧ Метод операторных рядов Чебышёва

МСРЛ Метод спектрального разложения Ланцоша

Таблица 2: Сокращения (русский алфавит).

Обозначение или сокращение Смысл или расшифровка

С Расширенная комплексная плоскость capl(-) Логарифмическая ёмкость card(-) Мощность множества

CJ Положительные константы (разные в разных главах)

Со(.) Выпуклая оболочка плоского множества diag{-, .,•} Диагональная матрица, составленная из заданных диагональных компонент diam(-) Диаметр множества в метрическом пространстве dist{-, •} Расстояние в метрическом пространстве ej единичный орт. Верхний индекс, если он есть, обозначает размерность пространства

FOU Метод полной ортогонализации

GMRES Обобщённый метод минимальных невязок

I Единичная матрица или тождественный оператор. Если есть индекс, то он обозначает размерность т-мерное подпространство Крылова

Mn(k) Кольцо матриц размера п х п над полем к

N, Q, R, С Множества соответственно натуральных, рациональных, вещественных и комплексных чисел

S Вещественная и мнимая части комплексного числа

Sp(-) Спектр матрицы или линейного оператора span{-, .,•} Подпространство, натянутое на заданные векторы

Supp(-) Носитель меры

Tk, Uk Многочлены Чебышёва первого и второго рода

Таблица 3: Обозначения и сокращения (латинский алфавит).

Краткое описание проблемы и основных результатов

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

Среди крыловских методов выделяются своей надёжностью (невозможностью преждевременных обрывов и переполнений, гарантией сходимости при определённых предположениях) классические методы Ланцо-ша и Арнольди, предназначенные соответственно для самосопряжённых и несамосопряжённых операторов. Однако существовавшая ранее теория сходимости методов Ланцоша и Арнольди была, на наш взгляд, в неудовлетворительном состоянии, так как содержала частные результаты или оценки не ошибок методов, а промежуточных величин.

В этой диссертации мы систематически излагаем построенную нами новую общую теорию сходимости известных и важнейших для практики методов Ланцоша и Арнольди, основываясь на полученных и опубликованных нами (частично с соавторами) в 1980-х — 2000-х годах результатах. Мы обсуждаем решение спектральных задач и вычисление матричных (операторных) функций, умноженных на вектор. По мере возможности мы изучаем поведение методов как в точной, так и в машинной арифметике.

На наш взгляд, в излагаемой теории наиболее существенны следующие достижения:

1. Систематическое использование поля значений (числового образа) матрицы (оператора) в оценках погрешности метода Арнольди как метода частичного вычисления спектра ограниченного оператора и как метода вычисления операторной функции, умноженной на вектор.

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

Стало: Оценки ошибок вычисления элементов спектра; достаточные условия сходимости (в операторном случае). Оценки погрешности вычисления операторных функций общего вида.

2. Получение улучшенных по порядку или не имевших аналогов оценок вычисления элементов спектра симметричных матриц методом Ланцоша и оценок вычисления умноженных на вектор функций от симметричных матриц методом Ланцоша в условиях машинной арифметики с помощью техники возмущённых чебышёвских рекурсий.

Было до нас: Оценки промежуточных величин. Сведение исследования процесса Ланцоша в машинной арифметике к исследованию процесса Ланцоша в идеальной арифметике с точностью до корня четвёртой степени из элементарной машинной ошибки округления.

Стало: Оценки погрешности вычисления методом Ланцоша в машинной арифметике матричных функций, квадратур и элементов спектра с точностью до умеренных кратных элементарной ошибки округления. Понимание того, что процесс Ланцоша в машинной арифметике может порождать на промежуточных шагах результаты, существенно различающиеся на разных компьютерах, но что это не мешает получению правильных финальных результатов.

3. Использование операторного спектра в оценках погрешности методов Ланцоша и Арнольди в точной арифметике.

Было до нас: Отдельные рассуждения об адаптации методов Ланцоша и Арнольди к спектру (зависимости результатов от спектра, позволяющей усилить оценки, основанные только на поле значений) и частные результаты о родственных методах.

Стало: Оценки погрешности методов Ланцоша и Арнольди в терминах операторного спектра. Открытие и объяснение явления сходимости на подпоследовательности шагов.

Результаты численных экспериментов подтверждают неулучшаемость многих наших оценок.

Основные результаты диссертации опубликованы в статьях [20, 21, 25, 26, 27, 28, 29, 90, 94, 114, 132], вышедших в журналах из списка ВАК. Другие результаты диссертации опубликованы в работах [91, 133]. Обзорные статьи, а также статьи с результатами, родственными вошедшим в диссертацию: [18, 89, 91, 94].

Ссылки на статьи автора по теме диссертации можно найти в книгах [11, предисловие], [110, § 7.5, 8.0], [113, гл. 4], [118, § 13.5-13.6], [142, § 3.11], [156, § 13.11], [182, гл. VI, § 28] и в статьях разных авторов.

19Автор полагает, что он был среди тех, кто в начале 1990-х годов стал активно внедрять использование поля значений, функции Грина и многочленов Фабера в исследование сходимости крыловских методов. Из «соратников» отметим М. Айерманна [97], О. Неванлинну [146] и JI. Н. Трефевена с соавторами [88] (список не претендует на полноту). Сейчас этими понятиями удивишь не очень многих вычислительных линейных алгебраистов: они (понятия) начали входить в учебники (см., например, [63, гл. 20]), а тогда их использование встречало заметное сопротивление.

Результаты диссертации докладывались на юбилейной лан-цошевской конференции [92, 128], конференции по вычислительной линейной алгебре [131], конференции 81АМ-САММ-2006 [134], на кафедре вычислительной математики Технического горного университета Фрайберга (Германия), в ВЦ РАН им. А. А. Дородницына, в ИПУ РАН им. В. А. Трапезникова, по нескольку раз обсуждались на семинарах в Институте вычислительной математики РАН и Институте прикладной математики РАН им. М. В. Келдыша (в том числе на семинаре им. К. И. Бабенко).

Часть I

Методы Ланцоша и спектрального разложения Ланцоша в точной арифметике: неадаптивные оценки

Глава 2

Метод спектрального разложения Ланцоша в точной арифметике

2.1 Введение

Пусть А — симметричная вещественная матрица размера п х п, / — функция, определённая на спектральном интервале А, <р £ Мп. Рассмотрим задачу вычисления вектора (1.7) и = /(А). Задачей такого типа является решение систем линейных уравнений (/(А) = А~1). В математической физике рассматриваемая задача возникает, например, при решении дифференциально-разностных аппроксимаций однородных уравнений в частных производных с коэффициентами, не зависящими от одной переменной (/(А) = ехр(—ЬА) — при решении параболических уравнений, /(А) = сов^Л1/2) — гиперболических уравнений, /(А) = ехр(—М1/2) — эллиптических уравнений1). Для вычисления и широко используются явные методы, в которых /(А) аппроксимируется многочленом (см., например, [13, 35, 38, 52, 143]). Рассмотрим общие свойства таких методов.

1 Связь матричных функций с дифференциальными уравнениями давно известна; см., например, [33].

Обозначим через А^, г = 1,2, . .,п, Лг+х > Л^, собственные значения А (которые в данном параграфе будем для простоты считать различными), а через — соответствующую ортонормированную систему собственных векторов. Если разложение <р имеет вид п v = фл,

2.1) то п

2—1

Возьмём в качестве приближённого значения и вектор ит = рт(А)(р^ где рт — многочлен степени не выше т — 1. Норма ошибки равна и - Um\\ — п

- Pm)(K)ViZi i=1

2.2) где г=1

6 — дельта-функция. Таким образом, задача сводится к полиномиальной аппроксимации / с обобщённым весом р (т. е. в узлах Лг- с весами <р\).

Основные вычислительные затраты в полиномиальных методах связаны с умножением А на векторы, поэтому главной характеристикой эффективности того или иного метода является число т, требуемое для получения заданной точности. Легко видеть, что при фиксированном т минимум в (2.2) достигается, когда рт — начальный кусок ряда Фурье / по многочленам, ортогональным с весом р. Процесс Ланцоша [138], [45, гл. 13], который первоначально предназначался для вычисления спектра А, позволяет в точной арифметике строить многочлены, ортогональные с весом рт, который в некотором смысле аппроксимирует р. Впервые в задачах типа (1.7) процесс Ланцоша был использован для решения линейных систем [154, 155, 171], затем в ряде работ он применялся для решения одномерных [147, 148] и многомерных [19] параболических и гиперболических уравнений. Общая схема применения процесса Ланцоша для решения задачи (1.7) была предложена в [18] ив другой форме — в [184].

Ниже получены оценки погрешности процесса Ланцоша при вычислении (1.7) в общем случае и для конкретных (наиболее употребительных) функций /. При доказательстве применяется разложение /(А) в матричный ряд Чебышёва, которое можно использовать как самостоятельный способ решения задачи (1.7). Применение рядов именно по многочленам Чебышёва обусловлено тем, что нужна оценка используемых многочленов не только на Эр (Л), но и в собственных значениях Н, а последние могут появляться всюду на спектральном интервале. Приведены также результаты численных экспериментов, показывающие сравнительную эффективность предлагаемых методов.

В этой главе изложение ведётся для точной арифметики. Возможность модификации утверждений из § 2-6 на случай машинной арифметики будет обоснована в гл. 4.

2.2 Метод операторных рядов Чебышёва

Будем считать, что матрица А не скалярная, т. е. Хп > Ах, и что = 1. В данном параграфе предполагается, что известны числа Ах и Ап, но результаты остаются справедливыми, если заменить Ах и Хп на их оценки снизу и сверху соответственно.

Положим

В = 5(ж) = ; лп — л\ Лп — Л1

Ап 4- Ах - (А„ - Ах)я где 1к — единичная матрица порядка к. Тогда ||В|| < 1 и функция д определена на отрезке [—1,1]. Рассмотрим ряд Чебышёва [46, гл. 3] функции оо д(х) = ^дкТк{х), (2.3) к=О где Тк — многочлены Чебышёва первого рода, ортогональные на интервале ] — 1,1[ с весом 1/\/1 — х2 и удовлетворяющие рекуррентному соотношению

То(я) = 1, Т1(х) = х, Тк+1(х) = 2хТк(х)-Тк-1(х), к> 1, (2.4) а дк — коэффициенты Чебышёва: j 1 -*-/ / а г ■'«г j ^ =- / -7===-ах. тг J \Д mm(2,k + l) j д(х)Тк{х) -1 ж 2

Поскольку и = /(А)(р = д(В)ср, то возьмём в качестве приближённого значения (1.7) вектор т—1 А;—О

Предложение 2.1 Пусть ряд (2.3) абсолютно сходится на отрезке [—1,1]. Тогда справедлива оценка погрешности оо и-Ут\\ < ]ГЫ <+00. (2.6) к=тп

Доказательство. Так как ТЦ1) = 1 для всех к и ряд (2.3) при х — 1 абсолютно сходится, то

00

22 ы < к—т

Учитывая неравенства ||1?|| < 1 и |ТЦж)| <1, — 1 < х < 1, получаем

00 оо п-ут\\<^Ы\\тк(в)\\ |М1<£Ы, к=т к=тп что и требовалось доказать. □

Предложенный метод вычисления (1.7), выражаемый формулой (2.5), назовём методом операторных рядов Чебышёва (МОРЧ). Несмотря на свою простоту, он даёт эффект, например, в эволюционных и некоторых других задачах, см. [179, 180]. Векторы Tk(B)ip могут быть вычислены рекуррентно, с помощью (2.4).

2.3 Метод спектрального разложения Ланцоша: оценка погрешности в терминах коэффициентов смещённого ряда Чебышёва вычисляемой функции

Напомним содержание т шагов процесса Ланцоша приближённого вычисления спектра оператора А. В подпространстве Крылова (1.1) /Ст(А, (р) строится базис ql,.,qm, получаемый в результате ортогонализации по Граму-Шмидту последовательности А(р,., Ат~1(р. Ортогонализацию можно провести с помощью трёхчленной рекурсии

Aqi = А-1дг-1 + ЗД + г = 1,2,., т - 1, (2.7) где poqo = 0, = у и Д > 0. Обозначим через Н трёхдиагональную матрицу

1 А \ А «2 02

Рт-1 ат через вг, вг = (^хг,., 5Шг)г, г = 1,2,.,т, — собственные значения и соответствующие нормированные собственные векторы Н\ е* есть г-й единичный вектор размерности т. Положим С} = ^ 1. Уг = ф^г

Тогда (9(,уг) — приближённые собственные пары оператора А, получаемые методом Ритца на /Ст(А, ф).

Предположим на время, что числа А; различны. Вследствие (2.1) тогда п и = У^гДАг)^. г=1

Ввиду ортогональности матрицы (5х. йт) ш т = 41 = = sliSi =

1 1=1 поэтому в качестве приближения ит к и естественно взять вектор т т г=1 г=1 что оправдывает определение (1.8) ит = Я/(Н)е1.

Положим ит = Я/(Н)е 1 вне зависимости от того, прост ли спектр А. Это определение корректно, поскольку спектр Н содержится в спектральном интервале А. Изложенный метод назовём методом спектрального разложения Ланцоша (МСРЛ).

Матричным выражением (2.7) является формула

АЯ-ЯН = РтЯт+А (2.8) частный случай (1.3) {АС} — (¿Н + Из неё следует, что если (Зт = 0 (это равенство должно стать справедливым по крайней мере на п-м шаге), то (9^ у{) — точные собственные пары А и ит = и.

Предложение 2.2 Пусть р £ ЩХ], с^р < т — 1. Тогда р(А)(р = Яр(Н)е 1.

Доказательство. Достаточно доказать с помощью конечной индукции, что

А'ч> = ЯН'е1, з<т-1. (2.9)

При ,7 = 0 равенство (2.9) принимает вид у? = и, очевидно, верно. Для перехода от у к ] + 1, используя (2.8) и (2.9) и учитывая, что, вследствие трёхдиагональности Н, е^Н^ех = 0 при 0 < ] < т — 2, получаем АА*<р = АЯ&ег = (ЯН + ртЯт+1еТт)Н3 ег = + (ЗтЧт+1еТгпН^1 = ЯН^е ь что и требовалось.2 □

Положим у 2 Ах Хп — Х\

2Отметим, что верен аналогичный результат (1.9) (/(А)<^ = Qf(H)e 1, / £ С[£], с^ / < тп — 1) для более общего метода — метода Арнольди. Фактически в доказательстве использована хессенберговость, а не трёхдиагональность Н. Отметим также, что этот простой результат был чуть позже опубликован Ю. Саадом в [164], но был известен последнему и раньше (частная переписка).

Поскольку, вследствие свойств аппроксимации Ритца,

Ai<6>¿<An, i = 1,2,., ш, (2.10) то \\V\\ < 1.

Теорема 2.1 Пусть ряд (2.3) абсолютно сходится на отрезке [—1,1]. Тогда справедливы равенство оо и - ит = ^ 9k[Tk(B)cp - QTk(V)ei} (2.11) k=m и оценка погрешности

00

11« - «mil < 2 Y; Ы < +оо. (2.12) k=m

Доказательство. Имеют место разложения в ряд оо u = J2gkTk(B)ip (2.13) k=О см. доказательство предложения 2.1) и

00 оо rn = Я^9кТк(У)ех = YJ9kQTk{V)el (2.14) к=0 к=О ряд сходится ввиду того, что ROOII < 1 и HQII = 1). Вычитая (2.14) из (2.13) и учитывая равенство Тк(В)ср = QTk(V)ei, к = 0,1,., m — 1 (предложение 2.2), получаем (2.11). Наконец, оо оо

11« - um\\ < Y, ы 01Т*(В)11 IMI + IIQII ||r*(V0||] < 2 Y, Ык=тп к=тп

Сравнивая (2.12) и (2.6), видим, что оценка погрешности метода спектрального разложения Ланцоша в два раза превышает оценку погрешности метода операторных рядов Чебышёва. Отметим, однако, следующее. На практике обычно погрешность метода операторных рядов

Чебышёва близка к правой части (2.6). Погрешность же метода спектрального разложения Ланцоша, как правило, существенно меньше правой части (2.12), что объясняется малостью величины в квадратной скобке выражения (2.11) при многих значениях к (при переходе от (2.11) и (2.12) квадратная скобка оценивалась тривиально). См. гл. 10 об адаптации к спектру.

12.1 Заключение

На наш взгляд, в диссертации систематически изложена теория сходимости методов Ланцоша и Арнольди и соответствующих методов вычисления произведения матричной функции и вектора. Доказаны многочисленные оценки погрешности решения частичной проблемы собственных значений и вычисления матричных функций, причём оценки для спектральных задач прямо или косвенно выведены из оценок для матричных функций — многочленов, выделяющих нужные собственные значения. По мере возможности освещено поведение методов как в точной, так и в машинной арифметике. Описана роль поля значений и операторного спектра; сформулирована концепция адаптации методов к (операторному) спектру; объяснено, что рассмотренные методы конкурентоспособны благодаря этой адаптации.

Не все приведённые результаты вполне удовлетворяют автора с эстетической точки зрения (см. список задач в § 12.2). Однако автор считает, что приведённых результатов достаточно для того, чтобы пользователи уверенно (не боясь неустойчивости) применяли методы Ланцоша и Арнольди.1

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

12.2 Задачи

Нам видятся следующие нерешённые задачи.

Глава 2, § 2.3, и гл. 10, § 10.4. Связанные вопросы: сколько нужно «заплатить» за переход в теореме 2.1 от рядов по многочленам Че-бышёва к рядам по ортонормированным многочленам, соответствующим паре (А, <р), и за ослабление условия теоремы 10.1, говорящего, что функция аналитична внутри «критической изолинии» функции Грина? Результаты § 10.2 и 10.3 убеждают в том, что «бесплатно» это сделать нельзя. При каких ограничениях на вычисляемую функцию можно получить оценки с нижним пределом и пределом максимума ошибки на двух последовательных шагах, аналогичные (10.1) и (10.9)?

Глава 4, теорема 4.1. Вопрос, аналогичный вопросу к теореме 2.1, но осложнённый машинной арифметикой.

Глава 5, теорема 5.1. Что будет, если взять экстремальный многочлен /, построенный в стиле [34], [36, § 11]? При этом многочлен уже не обязан принимать наибольшее (в пределах спектрального интервала) значение в искомом собственном значении. В случае точной арифметики ответ положителен (см. теорему 3.1).

Глава 6, предложение 6.1, теорема 6.1 и следствие 6.1. Вопрос про квадратуру Гаусса-Ланцоша, аналогичный вопросу к теореме 4.1.

В перечисленных выше задачах механический переход к другой системе базисных многочленов невозможен, так как требуется знать оценку значений многочленов в числах Ритца, которые могут «гулять» по спектральному интервалу, а также оценку соответствующих весов. На рабочей встрече с участием Э. Гринбаум, В. Л. Друскина и автора в середине 1990-х годов была высказана идея итерационно делать совместное уточнение указанных величин, начиная с неадаптивных к спектру многочленов (Чебышёва и т. п.) и в пределе заканчивая «наилучшими» многочленами. Насколько известно автору, эту идею никто не пытался реализовать.

Глава 9, теорема 9.1. Получить оценки погрешности вычисления геометрически кратного изолированного собственного значения и соответствующей системы корневых векторов.

Глава 10. Получить какие-то аналоги приведённых адаптивных результатов для фиксированной матрицы. Глава 11.

1. Теорема 11.1: доказательство локальной равномерности сходимости аппроксимант по 2.

2. Нахождение априорного условия представимости функции марковского типа, порождённой комплекснозначной мерой, в виде, удобном для применения метода Арнольди с ограниченным оператором (см. замечание 11.2). После решения этой задачи результаты предложений 11.3 и 11.4 автоматически распространятся на случай комплекснозначной меры.

3. Выяснение того, что можно сделать в случае рационального возмущения с невещественным или неположительным суммарным вычетом. Очевидно, рассматриваемая в данной статье квадратура ограничение вещественности и положительности снять не позволяет.

4. Объяснение видимой периодичности на рис. 11.3.

5. Оптимизация выбора матрицы в § 11.4.

12.3 Другие подходы и родственные результаты других авторов

12.3.1 Подход Э. Гринбаум к теории простого процесса Лан-цоша

Альтернативный (элегантный) способ получения оценок погрешности метода Ланцоша в условиях машинной арифметики предложила Э. Грин- 231 баум в работе [112]. Она для каждого фиксированного номера шага процесса т представила реальный процесс Ланцоша как процесс Ланцоша в точной арифметике, применённый к матрице большей размерности со спектром, лежащим на «малом» расстоянии от спектра А. К сожалению, расстояние приходится делать порядка корня четвёртой степени из элементарной ошибки округления, что приводит к худшим результатам, чем в [18, 21]. Однако некоторые вспомогательные результаты из [112] оказалось возможным использовать для получения более точных оценок.

12.3.2 Методы Ланцоша и Арнольд и и теория ортогональных многочленов

Теорема 8.1.11 из первого тома двухтомника [169, 170] утверждает, что если некая точка является изолированной точкой меры /л на единичной окружности Т, то подходящие корни соответствующих ортогональных многочленов сходятся к изолированной точке меры со скоростью геометрической прогрессии. Доказательство неконструктивно, так что знаменатель прогрессии не указан.

На самом деле метод Арнольди, применённый к оператору умножения на независимую переменную в пространстве Ь2,ц(Т) и к начальному вектору 1 (постоянной функции «единица»), генерирует соответствующие ортонормированные многочлены как векторы Арнольди, причём корни ортонормированных многочленов совпадают с числами Ритца. Так как оператор унитарен (и, значит, нормален) и так как круг строго выпукл, предложение 8.1 даёт конструктивную экспоненциальную оценку скорости сходимости. Асимптотически более точную (адаптивную) оценку можно извлечь из предложения 11.4.

Аналогично, теорема 1.1 из статьи [87] родственна лемме 10.3.

Сказанное не следует воспринимать как упрёк. Это скорее констатация того факта, что математики, работающие в разных областях, иногда фактически делают одно и то же, но не видят междисциплинарных связей из-за того, что говорят на разных математических «диалектах». В частности, не только математический анализ может быть полезен в вычислительной линейной алгебре, но и результаты из вычислительной линейной алгебры могут быть полезны аналитикам.

12.3.3 Псевдоспектр (спектральный портрет матрицы)

Понятие псевдоспектра, или, что почти то же самое, спектрального портрета матрицы, описано в книгах [10, 11, 182] (см. также библиографию в них). Псевдоспектр — прекрасное средство объяснения того факта, что нельзя прогнозировать поведение крыловских процессов с сильно ненормальными матрицами, основываясь только на спектре. Однако само по себе рассмотрение псевдоспектра не позволяет получать априорные оценки погрешности в удобных терминах, поскольку обычно отсутствует априорная информация о величине на изолиниях псевдоспектра значений многочленов с требуемой нормировкой. Кроме того, для получения оценок нужно ещё рассматривать псевдоспектры матриц-проекций Н.

Характерным в этом плане является § 28 "Arnoldi and related eigenvalue iterations" из книги [182], где приводится некая оценка погрешности вычисления собственного значения в методе Арнольди, сходимость правой части которой к нулю при т —¥ оо представляется неясной, а потом за дальнейшими результатами о сходимости чисел Ритца в методе Арнольди читатель отсылается к нашей статье [132].

12.3.4 Условные минимизационные задачи теории потенциала и адаптация к спектру

Недавно появились статьи (см., например, [76, 117, 136]), касающиеся приложений условных (с ограничениями) минимизационных задач теории потенциала к теории методов Ланцоша и Арнольди (только с нормальными матрицами). Полученные этой техникой линейно-алгебраические утверждения относятся не к фиксированному процессу Ланцоша или

Арнольди, а к последовательности процессов в Сп с растущим п,2 точнее, к тому, что в пределе при п —> оо получается на шаге тп процесса в Сп при условии тп/п —> Ь Е]0,1[. Это не то, что хотелось бы, но это какое-то приближение (отличное от нашего) к разгадке явления адаптации методов к спектру.

12.3.5 Методы с рестартами

Известно, что при применении метода Арнольди часто возникают трудности, связанные с возрастанием вычислительной работы (из-за реор-тогонализации) и возрастанием требований к памяти (из-за хранения векторов Арнольди) по мере увеличения номера шага. В простом процессе Ланцоша реортогонализация не делается, но приходится хранить векторы Ланцоша, если нужно вычислить сразу несколько матричных функций во всех или многих точках евклидова пространства (в геофизических приложениях это обычно не так). Стандартным лекарством от таких трудностей является рестарт процесса. По поводу развивающейся теории см., например, работы [67, 74, 80, 98, 181].

12.3.6 Итерации подпространства

Вместо нахождения индивидуальных собственных пар иногда (если индивидуальные собственные пары плохо обусловлены) удобнее итерационно найти приближение к натянутому на них инвариантному подпространству с помощью крыловского или блочного крыловского алгоритма (возможно, с рестартами). См., например, [63, гл. 9], [73]. Существуют также итерационные методы вычисления инвариантного подпространства, основанные на дихотомии спектра (см., например, [6, 43]).

2 Такие последовательности были введены ещё в [183].

12.3.7 Крыловские процессы с преобразованным оператором

Если вычисление матрично-векторной функции /(А)<р с помощью метода спектрального разложения Ланцоша или Арнольди отнимает много времени из-за плохой обусловленности или большой нормы матрицы А, то иногда имеет смысл применить тот же метод к матрице В = (/ + тЛ)-1, где т — параметр, подлежащий оптимизации. Естественно, что при этом соответственно меняется вычисляемая функция и что эффективность сильно зависит от «цены» умножения матрицы В на векторы (т. е. от времени факторизации матрицы I + тА или от времени итерационного решения системы линейных уравнений с этой матрицей. Примеры см. в [100, 144].

12.4 Благодарности

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

Автор благодарит Э. Гринбаум за совместное написание нескольких статей по теории метода Ланцоша, а также за полезные обсуждения личных статей автора.

Автор благодарит за полезные обсуждения диссертации в целом или части отражённых в ней статей и за библиографическую поддержку М. Айерманна, С. К. Асвадурова, Б. Бекерманна, А. В. Богатырёва,

Дж. Голуба, С. А. Горейнова, С. Н. Давыдычеву, Н. Л. Замарашкина,

X. Д. Икрамова, В. П. Ильина, Дж. Каллэм, А. В. Князева, В. И. Лебедева

О. В. Локуциевского, К. Любиха, Ю. М. Нечепуренко, Б. Парлетта, Л. Райхеля, А. Руэ, В. С. Рябенького, А. Скороходова, А. В. Собяни-на, В. 3. Соколинского, С. П. Суетина, К. Торреса-Вердина, Е. Е. Тыр-тышникова, М. Хохбрук, О. Эрнста, участников семинаров в Институте вычислительной математики РАН и Институте прикладной математики

РАН, участников конференций [92, 134].

Автор благодарит Дэвида Бейли за публикацию пакета [72] в Интернете.

Автор благодарит Центральную геофизическую экспедицию и Schlu-mberger-Doll Research за организационную помощь.

1. Абрамовиц М., Стиган И. М. Справочник по специальным функциям. М.: Наука, 1979. - 832 с.

2. Амосов А. А., Дубинский Ю. А., Копченова Н. В. Вычислительные методы для инженеров. М.: Высшая школа, 1994. 544 с.

3. Ахиезер Н. И. Классическая проблема моментов. М.: Физматлит, 1961. 311 с.

4. Бахвалов Н. С., Жидков Н. П., Кобельков Г. М. Численные методы. М.: Лаборатория базовых знаний, 2002. 630 с.

5. Бейтмен Г., Эрдейи А. Высшие трансцендентные функции. // Т. 2. Функции Бесселя, функции параболического цилиндра, ортогональные многочлены. М.: Наука, 1974. 296 с.

6. Булгаков А. Я., Годунов С. К. Круговая дихотомия матричного спектра // Сибирский матем. журнал. 1988. Т. 29. № 5. С. 59-70.

7. Воеводин В. В. Ошибки округления и устойчивость в прямых методах линейной алгебры. М.: Изд-во МГУ, 1969. 140 с.

8. Воеводин В. В., Кузнецов Ю. А. Матрицы и вычисления. М.: Наука, 1984. 319 с.

9. Гантмахер Ф. Р. Теория матриц. М.: Наука, 1988. 549 с.

10. Годунов С. К. Современные аспекты линейной алгебры. Новосибирск: Научная книга (ИДМИ), 1998. 389 с.

11. Годунов С. К. Лекции по современным аспектам линейной алгебры. Новосибирск: Научная книга (ИДМИ), 2002. 212 с.

12. Годунов С. К., Прокопов Г. П. Применение метода минимальных итераций для вычисления собственных значений эллиптических операторов // Ж. вычисл. матем. и матем. физ. 1970. Т. 10. № 5. С. 1180-1190.

13. Годунов С. К., Рябенький В. С. Разностные схемы. М.: Наука, 1973. 400 с.

14. Гончар А. А. О сходимости аппроксимаций Паде для некоторых классов мероморфных функций // Матем. сб. 1975. Т. 97 (139). № 4 (8). С. 607-629.

15. Гончар А. А., Суетин С. П. Об аппроксимациях Паде мероморфных функций марковского типа // Современные проблемы математики. 2004. Т. 5. С. 3-67.

16. Демидович Б. П., Марон И. А. Основы вычислительной математики. М.: Наука, 1970. 664 с.

17. Дзядык В. К. Введение в теорию равномерного приближения функций полиномами. М.: Наука, 1977. 512 с.

18. Друскин В. Л., Книжнерман Л. А. Использование операторных рядов по ортогональным многочленам при вычислении функций от самосопряжённых операторов и обоснование феномена Ланцоша. Деп. в ВИНИТИ. 1987. № 1535-В87. 47 с.

19. Друскин В. Л., Книжнерман Л. А. Спектральный дифференциально-разностный метод численного решения трёхмерных нестационарных задач электроразведки // Изв. АН СССР: Физика Земли. 1988. № 8. С. 63-74.

20. Друскин В. Л., Книжнерман Л. А. Два полиномиальных метода вычисления функций от симметричных матриц //Ж. вычисл. матем. и матем. физ. 1989. Т. 29. № 12. С. 1763-1775.

21. Друскин В. Л., Книжнерман Л. А. Оценки ошибок в простом процессе Ланцоша при вычислении функций от симметричных матриц и собственных значений // Ж. вычисл. матем. и матем. физ. 1991. Т. 31. № 7. С. 970-983.

22. Икрамов X. Д. Разрежённые матрицы // Итоги науки и техн.: Матем. анализ. Т. 20. М.: ВИНИТИ, 1982. С. 179-260.

23. Икрамов X. Д. Несимметричная проблема собственных значений. М.: Наука, 1991. 240 с.

24. Ильин В. П. Численный анализ. Часть 1. Новосибирск: Изд. ИВ-МиМГ, 2004. 335 с.

25. Книжнерман Л. А. Вычисление функций от несимметричных матриц с помощью метода Арнольди // Ж. вычисл. матем. и матем. физ. 1991. Т. 31. № 1. С. 5-16.

26. Книжнерман Л. А. Оценки погрешности метода Арнольди: случай нормальной матрицы // Ж. вычисл. матем. и матем. физ. 1992. Т. 32. № 9. С. 1347-1360.

27. Книжнерман Л. А. Качество аппроксимаций к хорошо отделённому собственному значению и расположение «чисел Ритца» в простом процессе Ланцоша // Ж. вычисл. матем. и матем. физ. 1995. Т. 35. № 10. С. 1459-1475.

28. Книжнерман Л. А. Простой процесс Ланцоша: оценки погрешности гауссовой квадратурной формулы и их приложения // Ж. вычисл. матем. и матем. физ. 1996. Т. 36. № 11. С. 5-19.

29. Книжнерман Л. А. Квадратура Гаусса-Арнольди для функции {(г1—А)~1<р, ср) и Паде-иодобная рациональная аппроксимация функций марковского типа // Матем. сб. 2008. Т. 199. № 2. С. 27-48.

30. Крейн С. Г. О классах корректности для некоторых граничных задач // Докл. АН СССР. 1957. Т. 114. № 6. С. 1162-1165.

31. Крылов А. Я. О численном решении уравнения, которым в технических вопросах определяются частоты малых колебаний материальных систем // Известия Академии наук СССР. VII серия. Отделение математических и естественных наук. 1931. № 4. С. 491-539.

32. Ланкастер П. Теория матриц. М.: Наука, 1978. 280 с.

33. Лаппо-Данилевский И. А. Применение функций от матриц к теории линейных систем обыкновенных дифференциальных уравнений. М.: Гос. изд-во технико-теор. л-ры, 1957. 454 с.

34. Лебедев В. И. Об итерационных методах решения операторных уравнений со спектром, лежащим на нескольких отрезках //Ж. вы-числ. матем. и матем. физ. 1969. Т. 9. № 6. С. 1247-1252.

35. Лебедев В. И. Явные разностные схемы с переменными шагами по времени для решения жёстких систем уравнений // Отдел вычислит. матем. АН СССР. 1987. Препринт № 177. 39 с.

36. Лебедев В. И. Функциональный анализ и вычислительная математика. М.: Физматлит, 2005. 296 с.

37. Ленг С. Алгебра. М.: Мир, 1968. 564 с.

38. Локуциевский В. О., Локуциевский О. В. Применение чебышёвских параметров для численного решения некоторых эволюционных задач // Ин-т прикл. матем. АН СССР. 1984. Препринт № 996. 30 с.

39. Марчу к Г. И. Методы вычислительной математики. М.: Наука, 1980. 536 с.40 4142