Содержание к диссертации
Введение
Глава 1. Интегральные уравнения трёхмерной стационарной задачи дифракции акустических волн 16
1.1. Классическая постановка трёхмерной стационарной задачи дифракции акустических волн 16
1.2. Интегральные формулировки исходной задачи 17
1.2.1. Система интегральных уравнений смешанного типа 17
1.2.2. Интегральные уравнения I рода с одной неизвестной плотностью 19
1.3. Обобщённые решения интегральных уравнений трёхмерной стационарной задачи дифракции акустических волн 20
1.3.1. Обобщённая постановка исходной задачи 20
1.3.2. Обобщённое решение смешанной системы интегральных уравнений 21
1.3.3. Обобщённые решения интегральных уравнений I рода с одной неизвестной плотностью 28
1.4. Вспомогательные утверждения 32
Глава 2. Интегральные уравнения краевых задач для уравнения Гельмгольца и их численное решение 41
2.1. Метод численного решения 41
2.1.1. Дискретизация интегральных уравнений 41
2.1.2. Аппроксимация интегральных операторов 43
2.1.3. О решении краевых задач на спектре интегральных операторов 52
2.2. Итерационные методы вариационного типа для решения СЛАУ с плотно заполненными матрицами 54
2.2.1. Основные понятия и обозначения 54
2.2.2. Критерии отбора итерационных методов 56
2.2.3. Алгоритмы итерационных методов и теоремы сходимости 57
2.2.4. Критерий остановки счёта 64
2.3. Численное решение краевых задач 65
2.3.1. Внутренние краевые задачи для уравнения Гельмгольца 65
2.3.2. Решение 1 краевой задачи на спектре 70
2.3.3. Численное исследование сходимости приближённых решений интегральных уравнений 1 и 2 краевых задач 72
2.4. Зависимость числа итераций от размерности СЛАУ 75
Глава 3. Численное моделирование дифракции акустических волн на трёхмерных включениях 81
3.1. Тестирование численного метода решения задач дифракции 81
3.1.1. Рассеяние плоской акустической волны на шаре 81
3.1.2. Проверка сходимости приближённых решений 85
3.1.3. Численное решение интегральных уравнений на спектре 85
3.2. Вычисление размерностей подпространств Крылова 93
3.3. Результаты компьютерного моделирования процессов дифракции 99
Заключение 109
Литература
- Система интегральных уравнений смешанного типа
- Обобщённое решение смешанной системы интегральных уравнений
- Дискретизация интегральных уравнений
- Рассеяние плоской акустической волны на шаре
Введение к работе
Математическое моделирование процессов распространения стационарных волн в средах с трёхмерными включениями играет важную роль в различных областях науки и техники и приводит к постановке достаточно сложных задач математической физики. Такие задачи принято называть задачами дифракции (трансмиссии) или задачами рассеяния. Они встречаются, например, в радиофизике, дефектоскопии, оптике, акустике океана и атмосферы, геофизике.
В данной работе рассматриваются, в основном, вопросы численного решения трёхмерных стационарных задач дифракции акустических волн. С математической точки зрения они заключаются в решении скалярных уравнений Гельмгольца, которые описывают процессы распространения акустических колебаний в трёхмерном пространстве и содержащемся в нём локальном включении. При этом искомые решения должны удовлетворять контактным условиям, заданным на границе включения, и условиям излучения на бесконечности. Кроме того, их поведение существенным образом зависит от отношения длин падающих волн к характерным размерам включений.
Аналитические решения задач дифракции могут быть найдены только в исключительных случаях, когда граница включения имеет достаточно простую геометрическую форму (например, сфера или эллипсоид). Для построения решений таких задач используются методы Винера-Хопфа, интегральных преобразований, степенных рядов [49, 50, 59, 62]. Важность аналитических методов заключается в том, что полученные с их помощью решения позволяют описывать исследуемые процессы с высокой степенью точности. Поэтому аналитические решения могут быть использованы в качестве тестовых при изучении задач дифракции с рассеивателями более сложной формы.
Ещё один подход к решению данного класса задач связан с изучением асимптотического поведения искомых решений при малых и больших значениях частоты исходного волнового поля. В случае низкочастотного рассеяния акустических волн для построения решений обычно используют методы теории
5 возмущений [40], а в случае высокочастотного рассеяния - лучевые методы,
метод параболического уравнения или метод эталонных задач [1, 64]. Однако
формальные преобразования, лежащие в основе этих методов, в большинстве
случаев не имеют строгого математического обоснования.
Следует отметить, что класс задач, охватываемый упомянутыми выше методами, остаётся весьма узким, поэтому основным методом исследования дифракционных процессов является прямое компьютерное моделирование.
Применение компьютера требует предварительного построения дискретного аналога (дискретной модели) исходной задачи, которое может быть выполнено различными способами. При сравнении получаемых различными способами дискретизации моделей предпочтение, очевидно, следует отдавать сочетающей в себе простоту и приемлемую точность описания исходной задачи. Такая модель может быть достаточно просто реализована в виде программы, предъявляя при этом минимальные требования к оперативной памяти и процессорному времени компьютера.
Дискретизацию граничных и гранично-контактных задач дифракции можно осуществить с помощью конечно-разностных и проекционно-сеточных методов [12, 17, 20, 38, 39, 41, 42, 45-47, 63, 71, 83, 89]. Как правило, их применение оправдано при решении задач, сформулированных в ограниченной области. Применение этих методов для построения дискретных аналогов внешних трёхмерных гранично-контактных задач связано с существенными трудностями, обусловленными тем, что расчётная область является неограниченной и трёхмерной. В этом случае неограниченную область приходится заменять её конечной подобластью, а условие излучения на бесконечности - краевым условием на внешней границе подобласти. Такая замена вносит в дискретную модель неустранимую погрешность, величина которой медленно убывает с увеличением подобласти, поскольку искомые решения исходных задач, как правило, медленно убывают на бесконечности.
Трудности усугубляются при рассмотрении задач дифракции с высокой частотой колебания, т.е. таких задач, в которых длина падающих волн мала по
сравнению с характерными размерами включения. В этом случае искомые решения сильно осциллируют, поэтому для достижения приемлемой точности при их аппроксимации необходимо уменьшать шаг сетки. Число узлов сетки в трёхмерных задачах изменяется как 0(h~3), где h - величина шага сетки, следовательно, уменьшение шага сетки приводит к быстрому росту размерностей дискретных моделей.
Перечисленные свойства трёхмерных задач дифракции приводят к тому, что алгоритмы их численного решения, полученные конечно-разностными и проекционно-сеточными методами, предъявляют слишком высокие требования к ресурсам компьютера и поэтому являются малоэффективными для решения данного класса задач.
Эффективным способом решения указанных проблем является переход от дифференциальной постановки исходной задачи к эквивалентной ей интегральной постановке, который может быть осуществлён методами теории потенциала. При этом трёхмерная задача в неограниченной области сводится к двумерной задаче, сформулированной на замкнутой поверхности включения. К недостаткам такого подхода следует отнести сложную для теоретического анализа структуру получаемых граничных интегральных уравнений, что компенсируется возможностью построения на их основе эффективных численных алгоритмов.
Переход к интегральным постановкам исходной задачи дифракции может быть осуществлён различными способами. Основная идея заключается в том, что отражённое и проходящее волновые поля ищутся в виде интефалов типа потенциала с неизвестными функциями, заданными на фанице включения. Ядрами этих интефалов являются фундаментальные решения соответствующих дифференциальных уравнений или их производные, поэтому они автоматически удовлетворяют как самим дифференциальным уравнениям, так и условию излучения на бесконечности для отражённого волнового поля. При этом исходная задача сводится к задаче отыскания таких неизвестных функций (плотно-
7 стей), которые обеспечат для искомых волновых полей выполнение контактных
условий на границе включения,
В данной работе для сведения исходных дифференциальных задач к эквивалентным им интегральным уравнениям систематически применяется непрямой вариант метода интегральных уравнений [5-9, 11, 15, 16, 29, 32, 34, 37, 48, 51-56, 66, 69, 70, 99], В этом варианте неизвестные плотности представляют собой вспомогательные источники волнового поля, распределенные по границе включения, тогда как в прямом варианте метода в качестве неизвестных функций выбираются граничные значения искомых волновых полей и их производные [32, 35, 98], Оба метода позволяют сводить исходную задачу к различным эквивалентным ей системам двух интегральных уравнений с двумя неизвестными плотностями. Важным достоинством непрямого метода интегральных уравнений является то, что он позволяет уменьшить число неизвестных функций и сформулировать исходную задачу в виде одного интегрального уравнения с одной неизвестной плотностью. Такие уравнения удобны для численного исследования и полученные в результате их дискретизации задачи предъявляют меньше требований к ресурсам компьютера по сравнению с другими эквивалентными формулировками.
Впервые подобный подход был предложен в публикациях [95, 96]. Наиболее законченные результаты его применения к исследованию стационарных задач дифракции в классических постановках изложены в работах [55, 88, 94].
В данной работе исходная задача рассматривается в обобщённой постановке. Для неё получена эквивалентная система двух интегральных уравнений с двумя неизвестными плотностями, а также пара различных слабо сингулярных интегральных уравнений I рода с одной неизвестной плотностью. Для этих уравнений проведено исследование условий эквивалентности исходной задаче и корректной разрешимости в классе обобщённых функций.
В отличие от интегральных уравнений II рода, методы численного решения которых хорошо разработаны [3, 12, 21, 42, 68, 81], численные методы решения интегральных уравнений I рода оставались в течение долгого времени
8 изученными весьма слабо. В общем случае задачи отыскания их решений являются некорректно поставленными. Однако наличие слабой особенности в ядрах соответствующих интегральных операторов позволяет использовать для приближённого решения таких уравнений численные методы, которые не содержат в явном виде процедуру регуляризации [3, 18, 36, 58, 61]. Это явление получило название «саморегуляризации» и в одномерном случае было исследовано в работах [5, 8, 9,13, 14, 67, 76, 77, 90, 103].
В данной работе развиты идеи, изложенные в работах [52-55]. Для приближённого решения исследуемых интегральных уравнений I рода используются алгоритмы, в которых неизвестная плотность отыскивается в виде линейной комбинации гладких финитных функций, образующих разбиение единицы на поверхности включения. При дискретизации интегрального уравнения поверхностные интегралы приближаются выражениями, содержащими интегралы по пространству R , которые затем вычисляются аналитически. Это позволяет рассчитывать коэффициенты систем линейных алгебраических уравнений, аппроксимирующих соответствующие интегральные уравнения, по весьма простым формулам. Такой подход, в отличие от наиболее популярного в настоящее время метода граничных элементов [2, 72, 84], не требует предварительной триангуляции поверхности, одинаково просто реализуется как на регулярных, так и на нерегулярных сетках, и позволяет обойтись без трудоёмкого приближённого вычисления кратных поверхностных интегралов. К его недостаткам следует отнести сложность теоретического обоснования. К моменту написания данной работы оно проведено только для одного вида эллиптических уравнений - уравнения Лапласа. Именно, в работе [55] показано существование и единственность приближённого решения интегрального уравнения I рода, которое эквивалентно внутренней и внешней задачам Дирихле для уравнения Лапласа, в классе обобщённых функций и получены оценки скорости убывания невязки и сходимости приближённого решения к точному.
Описание некоторых других подходов к численному решению интегральных уравнений, эквивалентных дифференциальным уравнениям эллиптического типа, имеется, например, в работах [6,7, 51, 69, 70, 74, 76, 78, 87,93].
Решение систем линейных алгебраических уравнений (СЛАУ), полученных в результате дискретизации интегральных уравнений I рода, можно находить различными методами. Поскольку основные матрицы таких систем являются плотно заполненными, количество операций, требуемых для нахождения решения прямыми методами, оценивается как 0{п ), где размерность матрицы я=10 -НО . Итерационные методы требуют 0{п ) операций на каждой итерации, поэтому, если число итераций невелико по сравнению с размерностью системы, использование итерационных методов позволяет существенно сократить время решения исследуемых задач.
Применяемые в данной работе методы отыскания приближённых решений СЛАУ являются итерационными методами вариационного типа, т.е. задача отыскания решения системы сводится к эквивалентной ей задаче минимизации некоторого функционала на каждом шаге итерационного процесса [73, 85, 102, 106]. Такой подход обладает следующим важным достоинством: он позволяет строить быстро сходящиеся итерационные процедуры и при этом не требует никакой априорной спектральной информации, за исключением сведений о невырожденности основной матрицы СЛАУ. С учётом того, что задача нахождения границ спектра квадратной матрицы общего вида является, вообще говоря, весьма сложной задачей, указанное достоинство оказывается очень существенным. Некоторые результаты применения одного из исследуемых методов (GMRES) для численного решения интегральных уравнений изложены в работах [75, 86].
После того, как приближённое решение интегрального уравнения найдено, искомое приближённое решение исходной задачи дифракции может быть одинаково просто и точно вычислено как в ближней, так и в дальней зоне.
Целью работы является разработка и реализация нового эффективного подхода к моделированию процессов дифракции стационарных акустических
10 волн в однородных средах с трёхмерными включениями. Исходя из поставленной цели, в работе рассматриваются следующие задачи:
получение и теоретическое исследование слабо сингулярных интегральных уравнений I рода с одной неизвестной функцией, эквивалентных задачам дифракции в обобщённых постановках;
разработка устойчивых численных алгоритмов с "саморегуляризацией" для решения полученных интегральных уравнений и их реализация в виде комплекса компьютерных программ;
математическое моделирование процессов дифракции стационарных акустических волн на трёхмерных включениях с использованием разработанного комплекса программ.
Методика исследований. Представленные в диссертации результаты теоретических исследований и вычислительных экспериментов получены с привлечением методов теории потенциала, дифференциальных и интегральных уравнений, теории обобщённых функций, функционального анализа и вычислительной математики.
Научная новизна работы состоит в следующем:
получены и исследованы новые обобщённые интегральные постановки стационарных задач дифракции акустических волн на трёхмерных включениях;
разработаны и реализованы в виде комплекса программ алгоритмы численного решения слабо сингулярных интегральных уравнений I рода с одной неизвестной функцией, эквивалентных задачам дифракции в обобщённых постановках;
разработан и численно реализован метод решения интегральных уравнений на спектре, т.е. в ситуациях, когда нарушаются условия эквивалентности дифференциальной и интегральной постановок исходных задач;
исследованы возможности применения итерационных методов вариационного типа для численного решения задач дифракции в интегральных постановках.
Теоретическая и практическая ценность. В диссертации получены и иссле-
дованы интегральные уравнения, позволяющие создавать эффективные численные алгоритмы решения трёхмерных задач дифракции.
Разработаны и реализованы в виде комплекса программ алгоритмы численного решения интефальных уравнений трёхмерных фанично-контактных задач акустики.
Проведены численные эксперименты, по результатам которых сделан вывод об эффективности итерационных методов вариационного типа для численного решения СЛАУ с плотно заполненными матрицами, аппроксимирующими интефальные уравнения трёхмерных задач дифракции.
Применяемые в данной работе методы могут быть использованы для численного исследования фаничных и контактных задач электродинамики и упругости путём их сведения к эквивалентным интефальным или интефо-дифференциальным уравнениям.
Публикации. Основные результаты диссертации опубликованы в 5 печатных работах.
Структура диссертации. Диссертация состоит из введения, трёх глав и заключения. Работа написана на 118 страницах и содержит 8 таблиц, 48 рисунков и список литературы из 106 наименований.
Перейдём к формулировке основных результатов диссертации.
Система интегральных уравнений смешанного типа
Аналитические решения задач дифракции могут быть найдены только в исключительных случаях, когда граница включения имеет достаточно простую геометрическую форму (например, сфера или эллипсоид). Для построения решений таких задач используются методы Винера-Хопфа, интегральных преобразований, степенных рядов [49, 50, 59, 62]. Важность аналитических методов заключается в том, что полученные с их помощью решения позволяют описывать исследуемые процессы с высокой степенью точности. Поэтому аналитические решения могут быть использованы в качестве тестовых при изучении задач дифракции с рассеивателями более сложной формы.
Ещё один подход к решению данного класса задач связан с изучением асимптотического поведения искомых решений при малых и больших значениях частоты исходного волнового поля. В случае низкочастотного рассеяния акустических волн для построения решений обычно используют методы теории возмущений [40], а в случае высокочастотного рассеяния - лучевые методы, метод параболического уравнения или метод эталонных задач [1, 64]. Однако формальные преобразования, лежащие в основе этих методов, в большинстве случаев не имеют строгого математического обоснования. Следует отметить, что класс задач, охватываемый упомянутыми выше методами, остаётся весьма узким, поэтому основным методом исследования дифракционных процессов является прямое компьютерное моделирование.
Применение компьютера требует предварительного построения дискретного аналога (дискретной модели) исходной задачи, которое может быть выполнено различными способами. При сравнении получаемых различными способами дискретизации моделей предпочтение, очевидно, следует отдавать сочетающей в себе простоту и приемлемую точность описания исходной задачи. Такая модель может быть достаточно просто реализована в виде программы, предъявляя при этом минимальные требования к оперативной памяти и процессорному времени компьютера.
Дискретизацию граничных и гранично-контактных задач дифракции можно осуществить с помощью конечно-разностных и проекционно-сеточных методов [12, 17, 20, 38, 39, 41, 42, 45-47, 63, 71, 83, 89]. Как правило, их применение оправдано при решении задач, сформулированных в ограниченной области. Применение этих методов для построения дискретных аналогов внешних трёхмерных гранично-контактных задач связано с существенными трудностями, обусловленными тем, что расчётная область является неограниченной и трёхмерной. В этом случае неограниченную область приходится заменять её конечной подобластью, а условие излучения на бесконечности - краевым условием на внешней границе подобласти. Такая замена вносит в дискретную модель неустранимую погрешность, величина которой медленно убывает с увеличением подобласти, поскольку искомые решения исходных задач, как правило, медленно убывают на бесконечности.
Трудности усугубляются при рассмотрении задач дифракции с высокой частотой колебания, т.е. таких задач, в которых длина падающих волн мала по сравнению с характерными размерами включения. В этом случае искомые решения сильно осциллируют, поэтому для достижения приемлемой точности при их аппроксимации необходимо уменьшать шаг сетки. Число узлов сетки в трёхмерных задачах изменяется как 0(h 3), где h - величина шага сетки, следовательно, уменьшение шага сетки приводит к быстрому росту размерностей дискретных моделей.
Перечисленные свойства трёхмерных задач дифракции приводят к тому, что алгоритмы их численного решения, полученные конечно-разностными и проекционно-сеточными методами, предъявляют слишком высокие требования к ресурсам компьютера и поэтому являются малоэффективными для решения данного класса задач.
Эффективным способом решения указанных проблем является переход от дифференциальной постановки исходной задачи к эквивалентной ей интегральной постановке, который может быть осуществлён методами теории потенциала. При этом трёхмерная задача в неограниченной области сводится к двумерной задаче, сформулированной на замкнутой поверхности включения. К недостаткам такого подхода следует отнести сложную для теоретического анализа структуру получаемых граничных интегральных уравнений, что компенсируется возможностью построения на их основе эффективных численных алгоритмов.
Переход к интегральным постановкам исходной задачи дифракции может быть осуществлён различными способами. Основная идея заключается в том, что отражённое и проходящее волновые поля ищутся в виде интефалов типа потенциала с неизвестными функциями, заданными на фанице включения. Ядрами этих интефалов являются фундаментальные решения соответствующих дифференциальных уравнений или их производные, поэтому они автоматически удовлетворяют как самим дифференциальным уравнениям, так и условию излучения на бесконечности для отражённого волнового поля. При этом исходная задача сводится к задаче отыскания таких неизвестных функций (плотно 7 стей), которые обеспечат для искомых волновых полей выполнение контактных условий на границе включения,
В данной работе для сведения исходных дифференциальных задач к эквивалентным им интегральным уравнениям систематически применяется непрямой вариант метода интегральных уравнений [5-9, 11, 15, 16, 29, 32, 34, 37, 48, 51-56, 66, 69, 70, 99], В этом варианте неизвестные плотности представляют собой вспомогательные источники волнового поля, распределенные по границе включения, тогда как в прямом варианте метода в качестве неизвестных функций выбираются граничные значения искомых волновых полей и их производные [32, 35, 98], Оба метода позволяют сводить исходную задачу к различным эквивалентным ей системам двух интегральных уравнений с двумя неизвестными плотностями. Важным достоинством непрямого метода интегральных уравнений является то, что он позволяет уменьшить число неизвестных функций и сформулировать исходную задачу в виде одного интегрального уравнения с одной неизвестной плотностью. Такие уравнения удобны для численного исследования и полученные в результате их дискретизации задачи предъявляют меньше требований к ресурсам компьютера по сравнению с другими эквивалентными формулировками.
Впервые подобный подход был предложен в публикациях [95, 96]. Наиболее законченные результаты его применения к исследованию стационарных задач дифракции в классических постановках изложены в работах [55, 88, 94].
Обобщённое решение смешанной системы интегральных уравнений
Пусть решение первой внешней краевой задачи для уравнения Гельмгольца ищется в виде интеграла типа потенциала простого слоя (1.16). Потенциал (1.16) удовлетворяет уравнению Гельмгольца в области С1е и условию излучения на бесконечности вида (1.3). Он будет решением внешней краевой задачи, если окажется возможным подобрать такую плотность q, которая обеспечит выполнение заданного граничного условия.
Полученное таким образом граничное интегральное уравнение эквива 53 лентно корректно разрешимой внешней краевой задаче при условии, что уе 0 или со не является собственной частотой задачи (1.7). Из леммы 1.12 следует, что при нарушении данного условия однородное граничное интегральное уравнение, в отличие от однородной краевой задачи, имеет нетривиальные решения. Поэтому неоднородное интефальное уравнение либо не имеет решения, либо имеет бесчисленное множество решений.
В этом случае, для нахождения с приемлемой точностью решения краевой задачи в точке спектра к интефального оператора или вблизи неё, будем поступать следующим образом: 1. Выберем одно или несколько таких волновых чисел, «близких» к волновому числу к исходной задачи, для каждого из которых алгоритм решения интегрального уравнения сохраняет устойчивость. 2. Решим интефальные уравнения с этими волновыми числами. 3. Вычислим решения соответствующих краевых задач в интересующих точках расчётной области. 4. Определим решение исходной краевой задачи как результат интерполяции полученных решений. Предложенный метод будем называть методом «интерполяции решения».
Пусть теперь рассматриваются сформулированные в виде фаничных интегральных уравнений внутренние краевые задачи для уравнения Гельмгольца. В этом случае спектр интефального оператора совпадает со спектром соответствующей краевой задачи, причём, если в точке спектра решение краевой задачи с неоднородными фаничными условиями существует, оно не является единственным.
Для такой задачи предложенный метод позволяет вычислять её частное решение, ортогональное ядру соответствующего дифференциального оператора. В 2.3.2 и 3.1.3 приведены результаты численных экспериментов, характеризующие возможности метода интерполяции решения для отыскания приближённых решений задач акустики на спектре интефальных операторов. 2.2. Итерационные методы вариационного типа для решения СЛАУ с плотно заполненными матрицами
Излагаемые ниже методы являются итерационными методами вариационного типа, т.е. исходная задача отыскания решения СЛАУ сводится к эквивалентной ей задаче минимизации некоторого функционала на каждом шаге итерационного процесса. Такой подход обладает следующим важным достоинством: он позволяет строить быстро сходящиеся итерационные процедуры и при этом не требует никакой априорной спектральной информации, за исключением сведений о невырожденности основной матрицы СЛАУ. С учётом того, что задача нахождения спектра квадратной матрицы является, вообще говоря, значительно более сложной задачей, чем задача нахождения численного решения СЛАУ с той же матрицей, указанное достоинство является весьма существенным. Рассмотренные итерационные методы применяются в последующих разделах для численного решения трёхмерных стационарных граничных и гранично-контактных задач дифракции в интегральных постановках.
Дискретизация интегральных уравнений
О выполнении или невыполнении критериев 2 и 3 можно судить по алгоритму итерационного метода, о выполнении критерия 1 - по результатам численных экспериментов. Что касается критерия 4, то существуют оценки, показывающие взаимосвязь невязок для некоторых методов. Эти оценки также будут использоваться при сравнении методов. Кроме того, о выполнении критерия 4 можно судить по результатам численных экспериментов.
Рассмотрев алгоритм метода минимальной невязки, можно сделать вывод, что метод удовлетворяет критериям 2 и 3. Скорость сходимости и число итераций, необходимых для нахождения решения методом минимальной невязки, существенно зависят от числа обусловленности решаемой системы. Поскольку рассматриваемые СЛАУ могут быть плохо обусловленными, критерий 1 может не выполняться.
Замечание 2.3. Метод минимальной невязки является одним первых методов вариациоиного типа. Раньше него был разработан только метод наискорейшего спуска, который позволяет решать СЛАУ с положительно опреде-лёнными матрицами и имеет одинаковую с методом минимальной невязки скорость сходимости. Мы не рассматриваем его, поскольку обобщение этого метода па случай систем с неположительными матрицами требует не одного, а двух умножений матрицы на вектор на каждой итерации, и потому заведомо уступает методу минимальной невязки.
Метод сопряжённых градиентов (Method of conjugate gradient, CG) был предложен в работе [82]. Данный метод является, по сути, прямым методом, поскольку позволяет получать точное решение СЛАУ с симметричной положительно определённой матрицей, при отсутствии ошибок округления, за п шагов, где п - порядок СЛАУ. В настоящее время метод используется как итерационный, поскольку во многих случаях для систем большого порядка число итераций, необходимых для нахождения приближённого решения СЛАУ с приемлемой точностью, значительно меньше п. Ниже приводится трёхчленная форма метода [80], которая эквивалентна его исходной форме.
Поскольку интегральные операторы трёхмерных стационарных задач дифракции акустических волн и их дискретные аналоги не обладают свойствами симметричности и положительной определённости, непосредственное применение метода сопряжённых градиентов для решения данного класса задач невозможно. Рассмотрим два обобщения метода сопряжённых градиентов на не симметричный невырожденный случаи. Первый подход изложен непосредственно в работе [82]. Он заключается в том, что обе части системы (2.22) умножаются слева на матрицу А . Полученная таким образом СЛАУ А Ах = А Ь (2.23) эквивалентна (2.22), её основная матрица А А симметрична и положительно определена. Таким образом, решение системы (2.22) можно найти методом сопряжённых градиентов.
Рассмотрев алгоритмы методов CGNR и CGNE, можно сделать вывод, что данные методы удовлетворяют критерию 3 и не удовлетворяют критерию 2, поскольку на одну итерацию в этих методах приходится два умножения матрицы на вектор. Скорость сходимости методов больше, чем у метода минимальной невязки, но следует учитывать, что число обусловленности основной матрицы системы в этих методах равно квадрату числа обусловленности исходной матрицы, поэтому для плохо обусловленных матриц скорость сходимости может существенно уменьшаться.
О выполнении или невыполнении критерия 4 будем судить по результатам численных экспериментов. Замечание 2.4. Рассматриваемые матрицы не являются нормальными, поэтому в системах (2.23) и (2.24) необходимо соблюдать порядок умножения матриц. Обобщённый метод минимальной невязки Обобщённый метод минимальной невязки (The Generalized Minimal Residual Method, GMRES) был предложен в работе [101]. Данный метод предназначен для решения СЛАУ с произвольной невырожденной матрицей, причём, как и в методе сопряжённых градиентов, при отсутствии ошибок округления, решение отыскивается не более чем за п шагов, где п - порядок СЛАУ.
Существует вариант GMRES, называемый перезапускаемый GMRES или GMRES(w), в котором указанный недостаток отсутствует. В этом варианте после каждых т итераций вычисляется промежуточное ре 63 шепне хт, которое затем принимается за повое начальное приближение х0 и процесс повторяется. У перезапускаемого GMRES существуют свои недостатки. Во-первых, в отличие от обычного GMRES, итерационный процесс сходится не для всех невырожденных матриц. Во-вторых, скорость сходимости у перезапускаемого GMRES, вообще говоря, меньше, чем у обычного GMRES. Скорость сходимости этого метода существенно зависит от числа т, максимальное значение которого определяется объёмом доступной оперативной памяти.
Рассеяние плоской акустической волны на шаре
Можно прогнозировать, что с ростом волновых чисел в задачах дифракции GMRES утратит свою эффективность, поскольку перестанет удовлетворять критерию 1 2.2.2. В этом случае необходимо применять GMRES с переобусловливанием, либо использовать его модификации, которые позволяют находить решение без использования всех базисных векторов. Кроме того, время решения задачи можно уменьшить за счёт выбора удачных начальных приближений и за счёт распараллеливания вычислений в многопроцессорных системах.
Общее время решения акустических задач дифракции можно оценить с помощью таблицы 3.3, в которой приведено время выполнения одной итерации в исследуемых итерационных методах для системы размерностью 8150. Таблица 3. итерационный метод CGNR CGNE GMRES затраты на одну итерацию, сек. 195 130 66 Замечание 3.2. Основная матрица СЛАУ размерности порядка 16000 уже не помещается целиком в оперативную память компьютера, что создаёт дополнительные сложности при реализации методов CGNR и CGNE, связанные с необходимостью дополнительно хранить на жёстком диске матрицу А . С этой точки зрения для решения СЛАУ большой размерности предпочтительно применять GMRES. Анализ алгоритмов итерационных методов и результатов вычислений показывает, что рассмотренные итерационные методы обладают рядом преимуществ по сравнению с прямыми методами при численном решении исследуемых интегральных уравнений, причём это преимущество значительно увеличивается с ростом размерности решаемых СЛАУ.
Во-первых, матрицы, аппроксимирующие исследуемые уравнения, являются плотно заполненными, следовательно, число операций, необходимых для решения СЛАУ прямыми методами оценивается как 0(п3). С учётом того, что количество итераций в рассматриваемых задачах, как правило, слабо зависит от размерности решаемой системы и на одну итерацию приходится 0(п2) операций, данный факт обуславливает преимущество итерационных методов при больших значениях п.
Во-вторых, для вычисления основной матрицы акустической задачи дифракции, необходимой при её решении прямыми методами, требуется допол-нительно ещё 0{п) операций (см. формулы (3.1)-(3.2)). Итерационные процедуры позволяют обходиться без вычисленной в явном виде основной матрицы и поэтому требуют лишь О(п ) операций на каждую итерацию.
Применение GMRES для решения интегральных уравнений трёхмерных стационарных краевых задач акустики позволяет сократить время вычислений (основная матрица порядка 6974) с 23,5 минут до 3,5-4 минут (здесь и далее - в зависимости от исходных данных), а непосредственно время решения СЛАУ-с 21 минуты до 1-1,5 минут.
При решении интегральных уравнений акустических задач дифракции с основной матрицей порядка 3998 применение GMRES даёт более ощутимый эффект и позволяет сократить время вычислений с 4,5 часов до 5-6 минут. Интересно, что при этом на непосредственно на решение СЛАУ методом Гаусса затрачивается 3 минуты, а остальное время уходит на формирование основной матрицы системы. Учитывая всё вышеизложенное, можно прогнозировать, что с ростом размерности СЛАУ преимущество итерационных методов решения будет только увеличиваться. Поэтому для данных задач переход от прямых методов решения к итерационным следует признать успешным, а итерационный метод GMRES рекомендовать для решения других задач математической физики, которые сводятся к граничным интегральным уравнениям.
Переходим к рассмотрению задач дифракции, которые отличаются от тестовой задачи из примера 3.1 видом включения, типом источника акустических волн и параметрами сред. В примерах 3.2 и 3.4 размерность СЛАУ выбиралась равной 16061, а в остальных примерах- 16033.
Приведённые графики распределения амплитуд колебаний акустического поля и модули диаграмм направленности позволяют делать выводы о расположении источника акустических колебаний относительно включения, размерах и форме отражающего объекта. Кроме того, видно, что максимум модуля добавочного поля ДФ находится вблизи включения, а ось JC2 является осью симметрии каждой из задач, кроме последней.
Численные эксперименты и приведённые примеры показывают, что используемый метод численного решения слабо сингулярных интегральных уравнений I рода обладает достаточно высокой точностью при решении интегральных уравнений трёхмерной стационарной задачи дифракции акустических волн. Он удобен для реализации на ЭВМ, не требует больших затрат машинного времени и может быть применён для решения других задач математической физики, которые сводятся к интегральным уравнениям по замкнутым поверхностям.
