Талтыкина (1235576), страница 4
Текст из файла (страница 4)
Элементы матрицыэтой линейной системы и использованная для дискретизации сетка являютсявходными данными для мозаично-скелетонного метода.2.2.1Построение дерева кластеровОбозначим через узел сетки с номером . Совокупность всех точек будем называть кластером нулевого уровня или корнем дерева. На каждомшаге кластер разбивается на непересекающиеся подкластеры, пока уровеньдерева не достигнет максимального значения.В качестве примера погрузим область, на границе которой построена сетка, в куб [21].
Все точки сетки, попавшие в куб, образуют кластер нулевогоуровня. На следующем шаге каждую сторону куба поделим пополам и получим 8 подкубов. Таким образом, все точки, принадлежащие исходному кубу,распределятся по 8 подкубам, образовав 8 подкластеров (возможно, некоторые из них будут пустыми).
Будем повторять операцию разбиения с подку25бами, получая при этом новые подкластеры. Если – последний уровеньразбиения, то количество кластеров этого уровня равно 8 .Обозначим центр куба -го уровня символом (или ), а сам куб символом ( ) (или ( )). Для любого куба -го уровня с центром в точке определим области или («зоны») () = { () ∈ : || − ||∞ 6 }, () = { () ∈ : 2( + 1) > || − ||∞ > ( + 1) },(2.12) – длина стороны куба уровня , - мозаичное разбиение, параметр небольшое целое число. Зоны будем называть «ближними», а – «дальними».Лемма 1. Пусть ∈ N – последний уровень иерархического разбиения.Тогда ∀, ∈ либо ∈ ( ), либо ∈ ( ) для некоторого номера, 1 < 6 .Доказательство леммы приведено в [21].Теорема 6.
Пусть матрица порождается ядром (1.8) на сетках (2.1),принадлежащих ограниченному множеству ∈ R3 , и подчиненных условию1mes Ωmes Ω 6 (Ω) 6 2,mes mes ∀Ω ∈ ,где (Ω) – число узлов сетки, принадлежащих области Ω. Пусть = ,где – диаметр множества . Тогда ∀, 0 < < 0 , и ∀ ∈ N, > 0 ,˜ для которыхсуществуют мозаично-скелетонные аппроксимации ,(︂)︂21mr ˜ < 1 + loglog ,˜ < 2 .|| − ||2.2.2Построение списка блоковПо дереву кластеров строится иерархический набор блоков разного размера, каждый блок отвечает взаимодействию двух кластеров.
При этом блоки,26отвечающие взаимодействию достаточно удалённых друг от друга кластеровточек, образуют «дальнюю» зону и могут быть приближены матрицей малого ранга. Блоки, описывающие взаимодействие рядом лежащих кластеров,определяют «ближнюю» зону. Элементы блоков «ближней» зоны вычисляются по формулам (2.4) для задач 1 и 2 и по формулам (2.9) для задачи3.2.2.3Построение малоранговой аппроксимацииВ работе [21] доказано, что существуют малоранговые аппроксимацииматриц, порождаемых осцилляционными ядрами интегральных операторов.К таким ядрам относится фундаментальное решение уравнения Гельмгольца(1.8).
Поиск малоранговой аппроксимации означает, что для каждого блока˜ для которой«дальней» зоны ∈ C× необходимо найти такую матрицу ,˜ ≤ ||||,˜|| − ||где – погрешность аппроксимации. При этом ранг матрицы ˜ должен бытьзначительно меньшим, чем её размеры. Для хранения блока в малоранговом˜ = ( + ) элементов, тогда как полныйформате нужна память mem()блок занимал бы mem() = элементов.
Малоранговая аппроксимацияпозволяет снижать объём памяти, необходимый для хранения каждого блокаи всей матрицы в целом [18].Малоранговую аппроксимацию можно построить разными способами, например, на основе сингулярного разложения матрицы или аналогичныхему методах типа Ланцоша. Такой подход использует все элементы матрицы, а также большое количество арифметических операций. Вследствие этого,на практике он редко применяется.В работах [15, 22] для этого предлагается находить в матрице столбцыи строки, на пересечении которых располагается подматрица максимальногообъёма (максимального модуля детерминанта) среди всех подматриц размера. Данная задача является весьма сложной в условиях неполной информациио матрице, поэтому чаще всего выбор подходящих столбцов и строк реализуется с помощью различных эвристических алгоритмов. Одним из таких27алгоритмов является метод Гаусса с выбором ведущего элемента по некоторому шаблону.
В работе [42] в качестве такого шаблона выбирается крестиз строки и столбца матрицы. По данному кресту на каждом шаге строитсятекущее приближение. Такой алгоритм называется неполной крестовой аппроксимацией, он является заключительным этапом мозаично-скелетонногометода [18].Алгоритм 1. Приближение матрицы размера × матрицей () , являющейся суммой скелетонов.Шаг 1. Пусть – номер вычисляемого скелетона. На этом шаге полагаем = 1 и выбираем произвольный номер столбца матрицы.Шаг 2. В столбце с номером вычисляем все элементы матрицы ,далее из них вычитаем элементы всех полученных ранее скелетонов вэтих позициях. В полученном столбце находим максимальный помодулю элемент.
Пусть он расположен в строке с номером .Шаг 3. Вычисляем строку матрицы с номером и из неё также вычитаем элементы всех уже найденных скелетонов в соответствующих позициях. В полученной строке * находим наибольший помодулю элемент, причём элемент из столбца с номером повторновыбрать нельзя. Номер столбца, в котором находится этот максимальный элемент, обозначим через +1 .Шаг 4.
По кресту с центром ( , ) строим скелетон так, чтобыкоэффициенты матрицы()=∑︁ *=1 (2.13)совпадали с коэффициентами исходной матрицы на позициях вычисленных столбцов 1 , ..., и вычисленных строк 1 , ..., .Шаг 5. Проверяем критерий остановки. Если он выполняется, товычисления прекращаются. Если не выполняется, то устанавливаем28 := + 1 и повторяем алгоритм с шага 2.Аппроксимация (2.13) считается достаточно точной, если выполняетсянеравенство (критерий остановки)|| − () || ≤ ||() || ,где – относительная погрешность аппроксимации, ||·|| – норма Фробениуса[43].Для проверки этого критерия нужно вычислить элементов матрицы, что слишком трудоёмко. Объём вычислений можно сократить, используяоценку сверху величины || − () || , полученную в [18],|| − () || ≤ ( − )|| * || /| | = ( − )|| || || || /| |,где = min(, ).При этом критерий остановки принимает вид( − )|| || || || /| | ≤ ||() || .(2.14)Норму матрицы () будем вычислять по рекуррентной формуле.
Для еёвывода воспользуемся определением нормы Фробениуса ||||2 = (* ),∑︀* – следгде – матрица, эрмитово-сопряжённая к матрице , () ==1матрицы [43]. Используя формулу (2.13), линейность следа и его инвариантность относительно циклической перестановки матриц [43], имеем:||() ||2(︁ ∑︁( * )( * ) )︁= = * ,=1=∑︁,=1(︁ ( * )( * ) )︁ * ∑︁(* )(* )=. * (2.15),=1Последнее выражение можно переписать в виде||() ||2=||(−1) ||2+ 2Re−1(︁ ∑︁(* )( * ) )︁=1* 29|| ||2 || ||2+,| |2(2.16)||(1) || =||1 || ||1 ||,|1 1 | ≥ 2.При использовании рекуррентной формулы (2.16) на -м шаге требуетсявыполнить (( + )) операций умножения.
Если в правой части (2.15) воспользоваться оценкой сверху для скалярных произведений, можно получитьболее предпочтительную с вычислительной точки зрения формулу()(−1)|| || ≤ ||∑︁ || || || |||| || || ||≤.|| +| |||=1(2.17)Она требует всего ( + ) операций на каждом шаге, однако при её использовании критерий остановки становится не вполне корректным, т.к. в правуючасть неравенства (2.14) вместо нормы ||() || подставляется её оценка сверху.
Вычислительные эксперименты, проведённые в широком диапазоне волновых чисел для решаемых ниже задач, показали, что использование формулы (2.17) не слишком ухудшает оценку точности малоранговой аппроксимации по сравнению с формулой (2.16), но позволяет уменьшить затратывремени на работу быстрого метода..303Программная реализация и численные эксперименты3.1Программное обеспечениеПрограммы для численного решения краевых задач и задач дифракциинаписаны на языке Fortran 90. В качестве компилятора используется IntelFortran Compiler, подключаемые библиотеки Intel Math Kernel Library (IntelMKL), а также реализация стандарта OpenMP 3.0 для данного компилятора.
Программное обеспечение представляют собой консольные приложения,которые работают на кластере ВЦ ДВО РАН. Для тестирования программбыл задействован вычислительный узел с 4-мя шестиядерными процессорамиSix Core AMD Opteronтм 8431, с тактовой частотой 2.40 ГГц и оперативнойпамятью 96 ГБ.Программы состоят из 21 модуля, которые реализуют все этапы решениязадач. Мозаично-скелетонный метод реализован в трёх модулях, соответствующих этапам метода. Так как краевые задачи формулируются в виде одного интегрального уравнения, то в этом случае программа находит решениеодной СЛАУ с одной матрицей, аппроксимирующей соответствующий интегральный оператор.
Для решения СЛАУ используются GMRES и мозаичноскелетонный метод. В случае задач дифракции матрица СЛАУ составляетсяиз четырёх матриц, аппроксимирующих операторы () , * и для решения задачи с помощью уравнения (1.15) или операторы () , * и длярешения задачи с помощью уравнения (1.17). Далее для простоты мы будем рассматривать программу решения задач дифракции, где будет СЛАУ сматрицей , аппроксимирующая оператор (2.4).Первый этап, построение дерева кластеров, реализован в виде рекурсивного алгоритма деления кластера, содержащего все точки сетки, на напересекающиеся подкластеры.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
