Main menu

Численная линейная алгебра: вычисление собственных значений и векторов

Скрытые резонансы: проблема спектрального анализа

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

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

Степенной метод для поиска доминирующего значения

Аналитически собственные значения находятся путем решения так называемого характеристического уравнения (поиска корней полинома, степень которого равна размерности матрицы). Однако, как мы знаем из теоремы Абеля-Руффини, для матриц размером 5x5 и больше точных формул корней не существует. Поэтому в реальности собственные значения всегда ищут только итерационными численными методами.

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

Алгоритм QR-разложения для нахождения всего спектра

Степенной метод хорош, если нам нужно только одно максимальное значение. Но что делать, если инженеру нужен весь спектр частот (все собственные значения)? Для этого в 1961 году Джоном Фрэнсисом и Верой Кублановской был независимо разработан QR-алгоритм, который считается одним из 10 самых важных численных алгоритмов 20-го века.

Алгоритм основан на факторизации (разложении) матрицы. На каждом шаге исходная матрица A представляется в виде произведения ортогональной матрицы Q и верхнетреугольной матрицы R (A = QR). Затем эти матрицы умножаются в обратном порядке: строится новая матрица A1 = RQ. Удивительный математический факт состоит в том, что новая матрица A1 имеет точно такие же собственные значения, что и исходная матрица A (они подобны). Если повторять этот процесс итеративно (находя QR-разложение и перемножая в обратном порядке), то матрица постепенно стремится к диагональному (или блочно-диагональному) виду, где все ее собственные значения просто стоят на главной диагонали. С применением предварительного приведения к форме Хессенберга и стратегий сдвигов, QR-алгоритм работает феноменально быстро и стабильно, являясь сегодня ядром библиотек линейной алгебры, таких как LAPACK, встроенных в Python, MATLAB и R.

Подробнее

Метод Монте-Карло: стохастический подход к вычислениям и интегрирование

Случайность на службе строгой математики

Численные методы, которые мы обсуждали до сих пор, являются детерминированными: они следуют строгим формулам, и при одинаковых входных данных всегда дают абсолютно идентичный результат с точностью до бита. Однако в середине 20-го века, во время работы над Манхэттенским проектом, физики Станислав Улам, Джон фон Нейман и Николас Метрополис предложили радикально иной подход, названный кодовым словом «Монте-Карло» в честь знаменитого казино. Этот метод использует генерацию случайных чисел для решения сугубо детерминированных математических задач.

Суть метода Монте-Карло заключается в проведении огромного числа виртуальных стохастических экспериментов. Классический пример — вычисление площади фигуры сложной формы или числа Пи. Если вписать фигуру в квадрат с известной площадью, а затем «бросать» случайные точки в этот квадрат, то отношение числа точек, попавших внутрь фигуры, к общему числу брошенных точек будет стремиться к отношению их площадей. Закон больших чисел теории вероятностей гарантирует, что при стремлении числа испытаний к бесконечности мы получим точное математическое решение.

Проклятие размерности в задачах интегрирования

Где же метод Монте-Карло проявляет себя лучше всего? Его главная ниша — вычисление кратных интегралов высокой размерности. В финансовой математике (ценообразование опционов), статистической физике или квантовой механике часто возникают интегралы, зависящие от десятков, сотен и даже тысяч переменных.

Если мы попытаемся использовать классические детерминированные методы (например, метод Симпсона) для вычисления D-мерного интеграла, мы столкнемся с «проклятием размерности». Если на каждую ось мы поместим по 10 узлов сетки, то для 3-мерного объема потребуется 10^3 = 1000 узлов. Но для 100-мерного интеграла потребуется 10^100 вычислений функции — число, превышающее количество атомов во Вселенной. Ни один суперкомпьютер с этим не справится. Ошибка классических методов зависит от размерности пространства. Метод Монте-Карло же обладает уникальным свойством: скорость его сходимости (убывания ошибки) равна O(1/sqrt(N)), где N — количество испытаний, и эта скорость абсолютно не зависит от размерности пространства! Это делает его единственным рабочим инструментом для задач высокой размерности.

Алгоритмы уменьшения дисперсии и цепи Маркова (MCMC)

Главный недостаток базового метода Монте-Карло — его медленная сходимость. Из-за корня в формуле сходимости, чтобы уменьшить ошибку в 10 раз, нужно увеличить количество испытаний в 100 раз. Чтобы бороться с этим, математики разработали изощренные техники снижения дисперсии (variance reduction techniques).

Одним из мощнейших инструментов является метод выборки по значимости (importance sampling). Вместо того чтобы разбрасывать случайные точки равномерно по всему объему, алгоритм концентрирует выборку в тех областях, где подынтегральная функция вносит наибольший вклад в итоговый результат. Другим прорывом стали алгоритмы Монте-Карло по схеме марковских цепей (Markov Chain Monte Carlo, MCMC), в частности алгоритм Метрополиса-Гастингса. В этом случае случайные блуждания не являются независимыми: следующая точка выбирается на основе положения предыдущей, формируя интеллектуальный поиск в многомерном пространстве параметров. Алгоритмы MCMC произвели революцию в байесовском машинном обучении, искусственном интеллекте и расшифровке генома, доказав, что контролируемая случайность — один из самых мощных инструментов познания.

Подробнее

Решение жестких систем обыкновенных дифференциальных уравнений (ОДУ)

Феномен жесткости: когда классические методы дают сбой

При численном моделировании многих химических, биологических и электронных систем исследователи сталкиваются с загадочным и крайне неприятным вычислительным явлением, получившим название «жесткость» (stiffness). Представьте себе химический реактор, в котором протекают две параллельные реакции: одна завершается за миллисекунды (например, горение или взрыв), а вторая длится часами (медленное окисление или диффузия). Эти процессы описываются системой обыкновенных дифференциальных уравнений (ОДУ).

Математически система называется жесткой, если ее решения содержат компоненты с кардинально различающимися скоростями затухания. Если мы попытаемся решить такую систему классическим явным методом (например, популярным методом Эйлера или методом Рунге-Кутты 4-го порядка), мы столкнемся с катастрофой. Для того чтобы обеспечить вычислительную устойчивость и не дать быстрым компонентам «взорвать» решение, шаг интегрирования должен быть меньше характерного времени самого быстрого процесса. В результате, чтобы промоделировать медленный процесс длительностью в часы, компьютеру придется сделать триллионы крошечных миллисекундных шагов, что займет месяцы машинного времени и приведет к огромному накоплению ошибок округления.

Неявные методы и А-устойчивость

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

На первый взгляд это кажется усложнением: теперь на каждом шаге нужно запускать метод Ньютона и решать системы линейных уравнений. Однако неявные методы обладают потрясающим свойством — абсолютной устойчивостью (или А-устойчивостью). Для таких методов размер шага интегрирования ограничивается только требованиями к точности описания физического процесса, а не условиями математической устойчивости. Медленный процесс можно интегрировать огромными шагами, полностью игнорируя сверхбыстрые затухающие колебания, не боясь при этом «взрыва» алгоритма. Неявный метод Эйлера является простейшим представителем этого класса.

Формулы дифференцирования назад (BDF) и методы Гира

Простейшие неявные методы обладают низким порядком точности (первым или вторым). Для точного моделирования жестких систем требуются методы более высоких порядков. Американский математик Чарльз Гир в начале 1970-х годов разработал семейство жестко-устойчивых алгоритмов, известных как формулы дифференцирования назад (Backward Differentiation Formulas, BDF).

Вместо того чтобы аппроксимировать подынтегральную функцию правой части, методы BDF строят интерполяционный полином по нескольким последним вычисленным точкам самого решения, а затем дифференцируют этот полином и приравнивают его к значению правой части в новой, еще не вычисленной точке. Это позволяет создавать высокоточные многошаговые неявные схемы (вплоть до 6-го порядка). Современные программные библиотеки (например, LSODE или решатели ode15s в MATLAB) используют адаптивные алгоритмы на основе методов BDF. Они автоматически распознают жесткость системы, на лету меняя как размер шага, так и порядок метода, обеспечивая робастное и эффективное решение сложнейших многомасштабных задач химической кинетики и схемотехники.

Подробнее

ГДЗ по математике

1 класс 2 класс 3 класс
4 класс 5 класс 6 класс

Книги

Фильмы

Старый сайт

Соц. сети