Файл: Российской федерации федеральное государственное автономное образовательное учреждение высшего образования.pdf

192
Предложенная реализация позволяет обмениваться моделями в виде интерактивных скриптов, просто передавая ссылку на модель. Возможные сценарии использования пакета, поддерживающие режим совместного ре- дактирования скриптов или допускающие только просмотр скрипта от- дельными участниками исследования.
Литература
1. Gerasimenko T. E., Kurbatova N. V., Nadolin D. K., Nasedkin A. V., Nased- kina A. A., Oganesyan P. A., Skaliukh A. S., and Soloviev A. N., Homogeni- zation of Piezoelectric Composites with Internal Structure and Inhomogene- ous Polarization in ACELAN-COMPOS Finite Element Package, in: Wave
Dynamics, Mechanics and Physics of Microstructured Metamaterials, edited by M. A. Sumbatyan (Springer, Berlin, 2019), pp. 113–131.
2. Официальный сайт библиотеки Prism.js. URL: https://prismjs.com (дата обращения: 10.01.2022).
3. Официальный сайт библиотеки
CodeMirror.
URL: https:// codemirror.net/6/.
4. Официальный сайт библиотеки
Earcut.
URL: https://github. com/mapbox/earcut.
5. M. Held. FIST: Fast Industrial-Strength Triangulation of Polygons. Algorith- mica 30(4): 563-596, 2001.
6. Официальный сайт библиотеки d3-delaunay. URL: https://github. com/d3/d3-delaunay.

193
РЕШЕНИЕ ЧАСТИЧНОЙ ПРОБЛЕМЫ СОБСТВЕННЫХ
ЗНАЧЕНИЙ ДЛЯ СИММЕТРИЧНЫХ
НЕЗНАКООПРЕДЕЛЕННЫХ МАТРИЦ НА
ПОДПРОСТРАНСТВАХ КРЫЛОВА
Мартынова Т. С., Муратова Г. В., Оганесян П. А., Штейн О. О.
ФГАОУ ВО «Южный федеральный университет»,
Институт математики, механики и компьютерных наук
им. И. И. Воровича
E-mail: martynova@sfedu.ru, muratova@sfedu.ru, poganesyan@sfedu.ru
Исследование выполнено за счет гранта Российского научного фонда
№ 22-21-00318, https://rscf.ru/project/22-21-00318/ в Южном федеральном университете.
Проблема собственных значений состоит в нахождении всех (полная проблема) или нескольких (частичная проблема) чисел λ
(исследуется вещественный случай), при которых существует ненулевое решение си- стемы Ax = λx, где A
–произвольная квадратная матрица. В работе рассмотрены методы на подпространствах Крылова, включающие соответ- ствующие рациональные методы, которые позволяют численно решать большие спектральные задачи для симметричных незнакоопределенных матриц, в частности, симметричных квазиопределенных (SQD) матриц.
Основными методами численного решения задачи являются методы Лан- цоша. Метод Ланцоша [1] первоначально являлся методом трехдиагонали- зации вещественной симметричной матрицы посредством ортогонального подобного преобразования. Позже он стал применяться во многих задачах: для вычисления собственных пар эрмитовых матриц [2, 3], как средство решения систем линейных алгебраических уравнений (СЛАУ), для локали- зации спектра в итерационных методах решения линейных систем и при решении жестких систем обыкновенных дифференциальных уравне- ний [4]. На практике метод Ланцоша используется для решения как знако- определенных, так и незнакоопределенных СЛАУ. Но при обосновании использования метода в последнем случае имеются трудности, связанные с проблемой устойчивости [5].
Пусть x – произвольный вектор n-мерного евклидова пространства с естественным скалярным произведением. Последовательность векторов
{x, Ax, A
2
x, ..., A
k
x, ...} называется последовательностью Крылова, порож- даемой матрицей A и вектором x. Линейные подпространства вида
K
m
= span{x, Ax, … A
m–1
x},
(1) называют подпространствами Крылова [2–4].

194
Особая привлекательность метода Ланцоша при решении больших спектральных задач состоит в том, что для проведения вычислений нет необходимости хранить A полным квадратным массивом. Необходимо лишь иметь возможность формировать произведения матрицы A на векто- ры, для чего матрицу можно хранить в компактной форме. Способ пред- ставления А определяет наилучший способ вычисления таких произведе- ний, но никак иначе не влияет на организацию алгоритма. Т.е. в алгорит- мах мы рассматриваем операцию матрично-векторного умножения как
«черный ящик». Открытый вопрос – это оптимальная для данной задачи конкретная реализация такого черного ящика. Это относится не только к разреженным матрицам, но и к заполненным матрицам, имеющим регу- лярную структуру, например, тѐплицевым или ганкелевым матрицам.
Будем искать несколько крайних (для определенности, максимальных) собственных чисел матрицы A, которую приведем к трехдиагональной форме Tm путем ортогонализации крыловской последовательности (1)
[2, 3]. Размерность матрицы Tm берется значительно меньшей размерности матрицы A (m << n). В [2] в качестве типичных значений для n, m и p, (где p – число искомых собственных чисел), приводятся p=10, m=300, n=104.
Величины p и m являются важнейшими параметрами методов Ланцоша.
Следует отметить, что числа Ритца (собственные значения матрицы Tm), которые считают приближениями к собственным значениям матрицы A, не обязаны лежать в ее спектре, они находятся в "поле значений" (множестве значений отношения Рэлея). Поэтому поиск собственных чисел матрицы, точный спектр которой неизвестен, представляет значительные трудности.
Однако известно [2, 3], что первыми сходятся крайние (максимальные и минимальные) собственные значения, причем сходимость тем быстрее, чем лучше они отделены друг от друга и от средней части спектра.
Пусть q0, q1, ... , qm-1 – ортобазис в подпространстве Крылова Km, и
Qm – матрица, столбцами которой являются данные векторы. Тогда она связана с матрицей A соотношением Tm = QТmAQm. Качество получен- ных приближений измеряется с помощью невязок [3].
Пусть µ – число Ритца, v = Qmy – соответствующий вектор Ритца (y – собственный вектор матрицы Tm, соответствующий собственному числу
µ). Известна оценка для нормы невязки пары Ритца [2, 3]
‖
‖
‖
‖
, (2) где
– последняя (m-я) компонента собственного вектора
, а
– единственный ненулевой элемент последней строки матри- цы
. Если норма невязки (2), окажется малой для неко- торых пар Ритца, то эти пары являются хорошим приближением соответ- ствующих собственных пар матрицы A. Заметим, что для вычисления

195 нормы невязки не нужно вычислять матрично-векторное произведение, а достаточно найти произведение двух чисел в правой части уравнения (2).
Однако, длина невязки, соответствующей вектору v,не является характе- ристикой его качества как приближенного собственного вектора без до- полнительной информации о спектре A. Собственные числа матрицы T
m
находятся с помощью стандартного QR-алгоритма.
Обычные методы подпространства Крылова аппроксимируют компо- ненты решения, связанные с доминирующими собственными значениями матрицы А, лучше, чем другие компоненты [2, 3]. Поэтому для жѐстких задач, где собственные значения сильно различаются по величине, может потребоваться много шагов Крылова m, прежде чем малые собственные значения матрицы А будут хорошоаппроксимированы. Обычно в данной ситуации используют два подхода: методы перезапуска или рациональные методы на подпространствах Крылова [5, 8]. Среди последних популярны методы типа «Shift-Invert», которые строят подпространство Крылова для матрицы (A – σI)
–1
, обратной к «сдвинутой» матрице A – σI. Следует отме- тить, что матрица (A – σI)
–1
не вычисляется. Вместо этого СЛАУ с матри- цей A – σI решается каждый раз, когда матрицу (A – σI)
–1
нужно умножить на вектор.
Если незнакоопределенная матрица А является симметричной квазио- пределенной (SQD) матрицей, то она является строго факторизуемой (име- ет модифицированное разложение Холецкого) [8], т.е. для любой матрицы перестановок Р существует нижняя треугольная матрица L и диагональная матрица D, такие что РТАР= LDLТ. При такой факторизации исключаются операции извлечения квадратных корней при вычислении диагональных элементов нижней треугольной матрицы. Кроме того, в случае LDLТ- разложения разреженной матрицы нижняя треугольная матрица L также является разреженной. Причем матрица L сохраняет характерные особен- ности структуры исходной матрицы, хотя общее количество ненулевых элементов, как правило, увеличивается. Поэтому для решения СЛАУ с матрицей A – σI можно использовать прямой метод. Альтернативой этому может служить решение данной СЛАУ алгебраическим многосеточным методом (или другим итерационным методом, учитывающим структуру и свойства SQD-матриц).
Следует отметить, что для нахождения собственных векторов задачи необходимо хранить весь набор векторов Ланцоша. Это требует большого объема памяти при работе алгоритма. Причем этот объем нельзя опреде- лить заранее, поскольку неизвестно число итераций, обеспечивающее вы- числение собственных значений с требуемой точностью. Для решения этой проблемы используются перезапуски метода Ланцоша [7].

196
Численные эксперименты
Первоначально тестирование программ производилось для диаго- нальной матрицы 500×500 (Таблицы 1 и 2), далее для модельной матрицы того же размера, полученной в результате конечно-элементной аппрокси- мации статической задачи теории электроупругости, спектр которой неиз- вестен (Таблица 3). В таблицах приводятся данные, отвечающие пяти мак- симальным собственным значениям данных матриц.
Численные эксперименты подтвердили, что метод Ланцоша находит экстремальные собственные значения с высокой точностью при достаточ- ном для сходимости выборе размерности подпространства Крылова
(m=50). При уменьшении размерности сходимость ухудшается (m=35), и точность вычисления собственных значений резко падает. Метод Ланцоша на рациональных подпространствах Крылова «Shift-Invert» показал высо- кую эффективность и точность при определении собственных значений, ближайших к заданному числу σ.
Таблица 1. Метод Ланцоша. Тестовая диагональная матрица, n=500, m=50, p=5
N Точное СЧ
Вычисленное СЧ
Погрешность
(абсолютная)
Невязка
1 2.814172761275914e+00 2.814172761275914e+00 0 1.55573281e-14 2
2.702868092389970e+00 2.702868092389970e+00 4.440892098e-16 8.19808960e-12 3
2.599089388576309e+00 2.599089388576309e+00 4.440892098e-16 2.47915414e-08 4
2.590000000000000e+00 2.589999999999998e+00 2.220446049e-15 4.87729831e-08 5
2.567541525093273e+00 2.567541525093271e+00 1.776356839e-15 3.50873753e-08
Таблица 2. Метод Ланцоша. Тестовая диагональная матрица, n=500, m=35, p=5
N Точное СЧ
Вычисленное СЧ
Погрешность
(абсолютная)
Невязка
1 2.814172761275914e+00 2.814172761275913e+00 8.881784197e-16 3.14414853e-08 2
2.702868092389970e+00 2.702868092373335e+00 1.663469362e-11 4.13708078e-06 3
2.599089388576309e+00 2.599080184787335e+00 9.203788973e-06 2.82704233e-03 4
2.590000000000000e+00 2.589972584614005e+00 2.741538599e-05 4.84766292e-03 5
2.567541525093273e+00 2.567534283211727e+00 7.241881546e-06 2.43875822e-03
Таблица 3. Метод Ланцоша. Модельная матрица, n=500, p=5
N
Вычисленное СЧ
(m=50)
Невязка
(m=50)
Вычисленное СЧ
(m=35)
Невязка
(m=35)
1 3.738423753468215e+00 1.51129110e-12 3.73842375e+00 9.83288595e-08 2
3.354251142333786e+00 9.83902222e-07 3.35424966e+00 1.35214213e-03 3
3.332708602332159e+00 9.82133017e-06 3.33262140e+00 1.01025964e-02 4
3.281254520234815e+00 6.01194574e-05 3.28039535e+00 2.96895087e-02 5
3.118811719900981e+00 6.71429252e-04 3.11858987e+00 8.23161078e-03

197
Таблица 4. Метод Ланцоша «Shift-Invert (σ)». Модельная матрица, n=500, p=5
N
σ
IT Вычисленное СЧ (λ(σ))
Ланцош (m=50) (λ
lanczos
)
Abs (λ(σ) - λ
lanczos
)
1 3.70 11 3.738423753491805e+00 3.738423753468215e+00 2.359001882723533e-11 2
3.35 13 3.354251142334455e+00 3.354251142333786e+00 6.687983500341943e-13 3
3.32 12 3.332708602401163e+00 3.332708602332159e+00 6.900391369413228e-11 4
3.25 16 3.281254523089236e+00 3.281254520234815e+00 2.854421143894115e-09 5
3.10 15 3.118812439169689e+00 3.118811719900981e+00 7.192687081492011e-07
В Таблице 4 приведены расчеты поиска методом Ланцоша «Shift-
Invert» пяти собственных чисел, ближайших к заданному числу σ. Резуль- таты сравниваются с соответствующими собственными значениями, полу- ченными методом Ланцоша, где «IT» обозначает число итераций, т. е. чис- ло обращений матрицы A – σI. Соответствующая СЛАУ решается моди- фицированным методом Холецкого, причем число итераций зависит как от близости σ к искомому собственному значению, так и от наличия двух
(или более) чисел Ритца, близко расположенных между собой.
Литература
1. Lanczos C., An iteration method for the solution of the eigenvalue problem of linear differential and integral operators, J. Res. Nat. Bur. Stand. В, v. 45, p. 255–281, 1950.
2. Parlett B. N., The Symmetric Eigenvalue Problem, SIAM, Philadelphia, PA,
1998.
3. Деммель Дж., Вычислительная линейная алгебра, М.: Мир, 2001.
4. Saad Y., Iterative methods for Sparse Linear Systems, PWS Publishing Com- pany, 1995.
5. Paige C. C., Parlett B. N., van der Vorst H. A., Approximate solutions and ei- genvalue bounds from Krylov subspaces, Numer. Linear Algebra with Ap- plic., v. 2, p. 115–133, 1995.
6. Ruhe A., Rational Krylov sequence methods for eigenvalue computation,
Linear Algebra Appl. v. 58, p.391–405, 1984.
7. Saad Y., Numerical methods for large eigenvalue problems, SIAM, second edition, 2011.
8. Vanderbei R. J., Symmetric quasidefinite matrices, SIAM J. Optim., 5(1), p. 100–113, 1995.

198