CC BY
12Научный вестник НГТУ. - 2010. - № 3(40)
УДК 519.6
Моделирование трехмерных магнитотеллурических полей векторным методом конечных элементов и возможности распараллеливания на процессорах с общей памятью*
П.А. ДОМНИКОВ, М.Г. ПЕРСОВА, Ю.Г. СОЛОВЕЙЧИК, Д.В. ВАГИН
Рассмотрены вопросы распараллеливания задач конечно-элементного моделирования трехмерных магнитотеллурических полей в общей памяти на уровне итерационных решателей систем линейных алгебраических уравнений. Использована вариационная постановка для векторного МКЭ с выделением основной части поля. Дано сравнение скорости сходимости итерационных методов BCG, COCG, COCR со сглаживанием невязки и GMRES для решения конечно-элементных систем уравнений, полученных в результате конечно-элементной аппроксимации трехмерной задачи магнитотеллурического зондирования.
Ключевые слова: метод конечных элементов, параллельное программирование, электромагнитные поля ВВЕДЕНИЕ
Задачи моделирования трехмерных геоэлектромагнитных полей в неоднородных средах - одни из наиболее трудоемких. Ускорение расчетов в этих задачах, несомненно, актуально, поскольку время проведения работ по восстановлению геоэлектрического строения среды зависит от вычислительных затрат, связанных с решением прямых трехмерных задач.
Для решения трехмерных задач геоэлектромагнетизма широко используется метод конечных элементов (МКЭ) [2, 3]. Использование специальных конечно элементных постановок с выделением основной части поля [1, 3] позволяет добиться очень высокой вычислительной эффективности и при моделировании трехмерных магнитотеллурических полей.
Рост возможностей вычислительной техники за последнее десятилетие позволяет разрабатывать методы решения сложных задач совершенно нового качественного уровня. Новые модели микропроцессоров, как правило, являются многоядерными, что применительно к задачам конечно-элементного моделирования дает возможность существенно ускорять процедуры получения численных решений за счет распараллеливания итерационных методов решения конечно-элементных систем линейных алгебраических уравнений (СЛАУ), в том числе и в общей памяти.
В работе рассматриваются вопросы распараллеливания методов решения СЛАУ, полученных при конечно-элементном моделировании трехмерных магнитотеллурических полей. Основное внимание будет уделено распараллеливанию итерационных методов решения конечно-элементных СЛАУ в многоядерных вычислительных машинах с общей памятью.
1. МАТЕМАТИЧЕСКАЯ МОДЕЛЬ И ВАРИАЦИОННАЯ ПОСТАНОВКА
Для моделирования трехмерного магнитотеллурического поля используется технология выделения основной части поля в горизонтально-слоистой среде [1, 3]. Вектор-потенциал магнитного поля А представляется в виде суммы поля А0 в горизонтально-слоистой среде (нор-
* Получена 16 марта 2010 г.
Работа выполнена при поддержке ФЦП «Научные и научно-педагогические кадры инновационной России» на 2009-2013 годы.
мальное поле) и добавочного поля Л+ от трехмерных неоднородностей (аномальное поле):
Л = Л0 + Л +.
Нормальное поле описывается следующей одномерной краевой задачей:
4 ' + 7О0юЛ0(z) = 0, (1)
ц0 dz2 1 dA0( z)
ц0 dz
= 1 + 0 • i. (2)
z=0
В (1)-(2) ^0 - магнитная проницаемость вакуума; с0 - удельная электрическая проводимость вмещающей горизонтально-слоистой среды; ю = 2tcv - циклическая частота поля; v - частота поля; i - мнимая единица. Под комплекснозначной величиной A0 понимается
x - или y -компонента вектор-потенциала A0. Для решения задачи (1)-(2) используется МКЭ с кусочно-линейными базисными функциями [3].
Для описания аномальной части поля в трехмерной области Q с границей используется векторное уравнение:
— rotrotA+ + ioraA+ = (о - о0)Ё0 в Q, A+ = 0 на <<Q, (3)
где с = с( x, y, z) - функция, описывающая удельную проводимость среды с включенными в
нее трехмерными неоднородностями, Ё0 = — iraA0 - нормальная составляющая напряженности электрического поля (в (3) она считается известной величиной). Индукция магнитного поля B и напряженность электрического поля Ё определяются соотношениями: B = rot A,
Ё = —<A , где A = A + + A0. <t
Эквивалентная вариационная постановка в форме Галеркина для задачи, определяемой уравнением (3), имеет вид
J—rot A+,Re • rotydQ — J o®A+,Im • уdQ = —J (о0 — o)mA0,Im • уdQ,
Q ^ Q Q
J—rot A+,Im • rotydQ + J ooA+,Re • уdQ = J (о0 — о)юЯ0,Re • уdQ,
(4)
+,Im + , f 7+,Re IV 0
-roL A ..../_
LQ ^0 Q Q
где у e H0 (rot,Q) = |у e [L2 (Q)]3 rot у e [L2 (Q)]3, [у x n]<Q = 0 j, т. е. H0 (rot,Q) - пространство вещественнозначных векторных функций, квадратично интегрируемых вместе со своим ротором в Q и имеющих нулевые тангенциальные компоненты на границе области .
Конечно-элементная аппроксимация уравнений (4) выполнялась на несогласованной сетке из параллелепипедов с кусочно-линейными векторными базисными функциями [3].
2. ОПИСАНИЕ ЗАДАЧИ МОДЕЛИРОВАНИЯ ТРЕХМЕРНОГО МАГНИТОТЕЛЛУРИЧЕСКОГО ПОЛЯ
Вычислительные эксперименты по сравнению итерационных методов решения СЛАУ и их распараллеливанию будем проводить на следующей трехмерной задаче магнитотеллуриче-ского зондирования (рис. 1).
Расчетная область состоит из вмещающей трехслойной среды с первым слоем толщиной 1500 м и сопротивлением 600 Ом-м, вторым слоем толщиной 220 000 м и сопротивлением 1000 Ом-м и фундаментом с сопротивлением 500 Ом-м. Размеры и проводимость находящихся в этой среде трехмерных неоднородностей приведены в табл. 1.
Таблица 1
Размеры и проводимость трехмерных неоднородностей
№ Размеры объекта, м Проводимость
п/п по X по у по г объекта, См/м
1 -20000...15000 5000...9000 -5500...-3500 0.01
2 -50000...-20000 -50000...50000 -15000...-5500 0.05
3 -20000..-5000 -5000...5000 -5500...-3500 0.01
Построенная нерегулярная конечно-элементная сетка содержит 37 707 параллелепипедов, 42 461 узел и 129 273 ребра. Фрагмент сечения конечно-элементной сетки в плоскости XZ показан на рис. 1.
Конечно-элементная дискретизация вариационной постановки (4) приводит к разреженной СЛАУ
Лх = /,
(5)
где комплексная матрица Л симметричная: (Л = Лт), но не эрмитова (Л Ф Лт).
Для данной задачи в конечно-элементной СЛАУ 106 861 комплексных неизвестных, что соответствует 213 722 вещественным неизвестным.
Рис. 1. Фрагмент сечения сетки в плоскости у = 0
2
3. ИТЕРАЦИОННЫЕ МЕТОДЫ РЕШЕНИЯ КОНЕЧНО-ЭЛЕМЕНТНЫХ СЛАУ И ИХ РАСПАРАЛЛЕЛИВАНИЕ В ОБЩЕЙ ПАМЯТИ
Матрицу СЛАУ вида (5) можно интерпретировать как вещественнозначную несимметричную и состоящую из блоков вида
ai] =
rb] —С] л
VС'] b'] J
(6)
где bj, Cj e M, и применять для решения СЛАУ (5) методы, ориентированные на решение
систем с несимметричными разреженными матрицами, например GMRES или BCG [5].
За последние годы были разработаны методы, позволяющие строить короткие рекурсии в подпространствах Крылова, если матрица СЛАУ является комплексно-симметричной (необязательно эрмитовой). К таким методам относятся, например, модификации метода сопряженных градиентов - Conjugate Orthogonal Conjugate Gradient Method (COCG) [7] и сопряженных невязок - Conjugate A-Orthogonal Conjugate Residual Method (COCR) [6]. Формулы этих методов очень похожи на формулы методов сопряженных градиентов и сопряженных невязок для симметричных вещественнозначных матриц. Так, формулы метода COCG имеют следующий вид:
Вычислить r0 = f - Ax0, p0 = r0 , (7)
а, = -_ J—, (9)
] (Ap,,P, )
Для j = 0,1,2,... (8)
jj)
pj, Pj
Xj+1 = XJ + аJPJ , (10)
rj+1 = rj _ а JAPJ , (11)
P j = Jj2, (12)
J (rj, rj )
PJ+1 = rj+1 + PJPJ . (13)
Здесь A e Cnxn, Xj,bj, Pj, rj e С", аj ,вj e С Vj. Черта сверху означает комплексное сопряжение. Скалярное произведение определяется формулой: (x,y) = ^"=1 .
В методе COCR коэффициенты ау- и вв (9), (12) вычисляются по формулам: а j = (rj, Arj) / (Apj, Apj), в j = (rj+1, Arj+1) / (rj , Arj); вектор Apj пересчитывается рекур-рентно: Apj = Arj + в j _1 Apj _1.
Невязки методов COCG и COCR подвержены сильным осцилляциям, что затрудняет контроль сходимости методов и завершение процесса по исчерпанию числа итераций. Поэтому для данных методов целесообразно использовать процедуру сглаживания невязки (minimal residual smoothing - MRS), обобщение которой для комплекснозначных СЛАУ приведено в [4].
Суть сглаживания невязки заключается в следующем [8]. Пусть некоторый итерационный метод генерирует на j-й итерации приближенное решение Xj и соответствующий вектор
невязки rj . Вводятся два вспомогательных вектора y j и s j следующим образом:
y0 = x0, s0 = r0,
У] = (1 -Пз)Уз-1 + Пзхз, = (1 - Пз )5У-1 + П3Г3 для ] = 1,2,3,...
II 1|2
Коэффициент пз е К выбирается так, чтобы величина / — АуЛ минимизировалась на каждом шаге:
- ('з-1, гз - *з-1) + (гз - *з-1,-1) з 2(гз - 5з-1, гз - 5з-1) '
Решением системы на з'-й итерации считается вектор Уз с соответствующим вектором
невязки Sj (||<|Гз||).
Сравнение сходимости метода СОСG без сглаживания и со сглаживанием невязки для задачи с частотой V = 0.1468 Гц приведено на рис. 2.
1E+0.2 1E+0.1
s
I 1E+0.0
£В К
'§ 1E-0.1
¡г
hC
I 1E-0.2-
К
О
Is 1E-0.3
о
rt
% 1E-0.4
Я
1E-0.5
0 1000 2000 3000 4000 5000
Число итераций
Рис. 2. Сходимость метода COCG со сглаживанием невязки
Модификация итерационного метода с использованием сглаживания невязки обладает более плавной сходимостью и в некоторых случаях позволяет достигать требуемой точности за меньшее число итераций, чем в исходном методе. В рассмотренном примере метод COCG
уменьшил норму относительной невязки в 104 раз за 5429 итераций, а COCG(MRS) - за 4225 итераций.
В табл. 2 приведено сравнение скорости сходимости методов GMRES, BCG, COCG и COCR для задачи моделирования трехмерного магнитотеллурического поля в расчетной области, изображенной на рис. 1, при различных значениях частоты поля v . В таблице указано число матрично-векторных умножений (как самых затратных в вычислительном плане операций в решателях) и время решения СЛАУ на персональном компьютере c микропроцессором Intel Core 2 Quad Q6600 2.4 GHz. Метод GMRES имел глубину 50. Матрица СЛАУ хранилась в разреженном строчно-столбцовом формате. Процесс решения заканчивался при уменьшении нормы относительной невязки в 104 раз. Для всех методов, кроме GMRES, использовалась технология сглаживания невязки (метод GMRES уменьшает невязку на каждой итерации сам по себе). В качестве предобусловливателя использовалась главная диагональ матрицы, состоящая из комплексных чисел. Прочерк в таблице означает, что данным методом решить СЛАУ не удалось - было превышено 100 000 матрично-векторных умножений.
Таблица 2
Число матрично-векторных умножений (и) и время решения СЛАУ (¿), требуемые для решения задачи различными итерационными методами (последовательными версиями) в зависимости
от частоты поля
V , Гц 0МЕЕБ(50) ВСО(МЕ5) СОСО(МЕ5) СОСЩМЯБ)
п t п t п t п t
0.000584 - - - - 27239 17 мин 20 с 20217 13 мин 14 с
0.001166 - - - - 22342 14 мин 9 с 16331 10 мин 43 с
0.002326 - - - - 17553 11 мин 5 с 13038 8 мин 33 с
0.004642 - - - - 13569 8 мин 35 с 10196 6 мин 41 с
0.009261 - - - - 11655 7 мин 21 с 8133 5 мин 20 с
0.01848 - - - - 9238 5 мин 50 с 6383 4 мин 10 с
0.03687 - - - - 7264 4 мин 35 с 5112 3 мин 21 с
0.07356 - - - - 5254 3 мин 19 с 4105 2 мин 42 с
0.1468 - - - - 4238 2 мин 41 с 3238 2 мин 6 с
0.2929 - - - - 3178 2 мин 1 с 2563 1 мин 41 с
0.5843 - - - - 2367 1 мин 29 с 1810 1 мин 11 с
1.166 - - 74766 50 мин 1 с 1553 59 с 1237 49 с
2.326 - - 28814 19 мин 24 с 1023 39 с 926 36 с
4.642 50650 59 мин 31 с 13678 9 мин 12 с 721 28 с 602 24 с
9.261 13650 16 мин 12 с 9282 6 мин 16 с 471 18 с 450 18 с
36.87 1200 1 мин 25 с 2546 1 мин 43 с 230 9 с 207 8 с
146.8 200 14 с 748 31 с 131 5 с 81 4 с
316.2 100 7 с 310 12 с 87 4 с 58 3 с
По приведенным в табл. 2 данным можно заключить, что, во-первых, время решения данной задачи рассмотренными итерационными методами возрастает с уменьшением частоты V электромагнитного поля, во-вторых, методы GMRES и BCG очень чувствительны к понижению частоты электромагнитного поля (начиная с некоторой частоты число итераций резко возрастает и при дальнейшем понижении частоты данные методы перестают сходиться). В-третьих, СОСЯ требует наименьших вычислительных затрат среди рассмотренных методов.
Теперь перейдем к вопросам распараллеливания итерационных методов в системах с общей памятью. В данной работе для распараллеливания в общей памяти использовались директивы ОрепМР. В рассмотренных итерационных методах могут быть распараллелены следующие основные операции: скалярное произведение и линейное комбинирование векторов, диагональное предобусловливание и умножение разреженной матрицы на вектор [5]. Остановимся на вопросах распараллеливания матрично-векторного умножения, поскольку именно оно определяет основные вычислительные затраты решателя.
Рассмотрим следующие два формата хранения разреженных матриц в памяти компьютера [3]. В разреженном строчном формате (РСФ) все ненулевые элементы матрицы (включая элементы из верхнего и нижнего треугольника матрицы) хранятся последовательно в строчном формате в одном массиве. В разреженном строчно-столбцовом формате (РССФ) для симметричной матрицы внедиагональные элементы матрицы хранятся только для нижнего треугольника (перечисленные в строчном порядке). При таком способе хранения элементы верхнего треугольника оказываются расположенными в том же массиве, но в столбцовом порядке. Положение ненулевых элементов в матрице СЛАУ в рассматриваемой задаче моделирования трехмерного магнитотеллурического поля и распределение элементов матрицы по процессорам (ядрам) для различных форматов хранения показаны на рис. 3.
а б в
Рис. 3. Расположение ненулевых элементов в матрице СЛАУ (а), распределение элементов матрицы по четырем процессорным ядрам для разреженного строчно-столбцового формата (б) и для разреженного строчного формата (в)
Преимущество разреженного строчно-столбцового формата заключается в том, что требуется в два раза меньше памяти для хранения внедиагональных элементов матрицы (и соответственно в два раза меньше обращений к памяти), чем при использовании разреженного строчного формата. При хранении матрицы в разреженном строчно-столбцовом формате также упрощается умножение на транспонированную матрицу. С другой стороны, при распараллеливании операций матрично-векторного умножения с матрицей, хранящейся в разреженном строчно-столбцовом формате, могут возникать конфликты при одновременной записи в вектор-результат (эта же проблема возникает и при умножении вектора на транспонированную матрицу). Устранить такие конфликты можно следующим образом: каждая нить параллельной программы производит запись в «свой» вектор-результат, а затем сложение всех векторов-результатов в один.
При использовании разреженного строчно-столбцового формата могут возникать конфликты и при одновременном чтении разными нитями параллельной программы компонент вектора, на который умножается матрица (т. е. при чтении входного вектора). В отличие от конфликтов при записи в вектор-результат, конфликты при чтении не приводят к ошибкам, так как компоненты входного вектора остаются постоянными. Однако одновременное обращение разных нитей к одной области памяти может замедлить работу программы, поэтому при распараллеливании матрично-векторного умножения с матрицей, хранящейся в разреженном строчно-столбцовом формате, имеет смысл провести эксперимент, в котором будут сделаны копии входного вектора для каждой нити (рис. 4).
Одна из существенных проблем распараллеливания операции умножения разреженной матрицы на вектор - неравномерная заполненность матриц. При распределении строк разреженной матрицами между нитями параллельной программы процессоры (ядра) могут быть
Рис. 4. Иллюстрация к распараллеливанию матрично-векторного умножения для матрицы, хранящейся в разреженном строчном формате (вверху) и в разреженном строчно-столбцовом формате (внизу). Операции, указанные на рисунке в круглых скобках, можно
выполнять параллельно
загружены неравномерно, что может вызвать снижение общей производительности. Для устранения возможных простоев строки матрицы распределяются по процессорам (ядрам) таким образом, чтобы каждому из них досталось приблизительно равное количество элементов матрицы. Пример распределения элементов матрицы для задачи моделирования трехмерного маг-нитотеллурического поля при распараллеливании на 4 процессорах (ядрах) с равномерным и неравномерным (т. е. с выравниванием по числу элементов) распределением строк приведен в табл. 3.
Таблица 3
Распределение элементов матрицы при использовании различных форматов
хранения
Формат матрицы, распределение строк Номе] р нити
1 2 3 4
РССФ, неравномерное 890 279 888 304 888 340 883 654
РССФ, равномерное 806 855 865 361 866 783 1 011 578
РСФ, неравномерное 887 693 887 630 887 655 887 599
РСФ, равномерное 897 503 895 409 899 115 858 550
Таблица позволяет заключить, что при использовании разреженного строчно-столбцового формата с равномерным распределением строк возникает ощутимый дисбаланс распределения элементов между нитями параллельной программы: 1-я нить обрабатывает 806 855, а 4-я - 1 011 578 элементов матрицы. При неравномерном распределении строк количество элементов матрицы, обрабатываемых каждой нитью программы, становится почти одинаковым.
Сравнение производительности решателя COCR(MRS) при использовании различных форматов хранения разреженных матриц на примере решения СЛАУ, полученной из задачи конечно-элементного моделирования трехмерного магнитотеллурического поля на частоте 0.000584 Гц, приведено в табл. 4. Распределение строк по ядрам в этом эксперименте было сделано с выравниванием числа элементов, обрабатываемых каждым ядром.
Таблица 4
Сравнение производительности решателя COCR(MRS) при распараллеливании с использованием разных форматов хранения разреженной матрицы СЛАУ
Формат хранения матрицы Число ядер
1 2 3 4
время, с время, с ускорение время, с ускорение время, с ускорение
РССФ 13 м 14 с 8 м 37 с 1.536 6 м 37 с 2.000 6 м 14 с 2.123
РССФ с копированием входного вектора 13 м 14 с 8 м 22 с 1.582 7 м 3 с 1.877 7 м 16 с 1.864
РСФ 14 м 36 с 9 м 2 с 1.616 7 м 19 с 1.995 6 м 53 с 2.121
Как видно из результатов, приведенных в табл. 4, разреженный строчно-столбцовый формат выигрывает по времени решения СЛАУ и по ускорению за счет распараллеливания у разреженного строчного формата. К тому же разреженный строчно-столбцовый формат требует в два раза меньше памяти при хранении симметричных матриц.
Копирование входного вектора в подпрограмме матрично-векторного умножения привело к увеличению ускорения только при распараллеливании на два ядра. При дальнейшем увеличении числа ядер ускорение, наоборот, уменьшалось, очевидно, из-за накладных расходов на копирование входного вектора.
Таким образом, при распараллеливании в общей памяти предпочтение следует отдать разреженному строчно-столбцовому формату. Следует заметить, что данный вывод справедлив при использовании небольшого числа процессоров/ядер: 2, 3, 4.
Ниже приведем время решения конечно-элементных СЛАУ, полученных из конечно-элементной дискретизации трехмерной задачи магнитотеллурического зондирования на разных частотах параллельной версией метода COCR(MRS) как наиболее быстрого метода решения данной задачи. Матрица СЛАУ хранилась в разреженном строчно-столбцовом формате, распределение элементов матрицы между ядрами оптимальное. При распараллеливании на два ядра в подпрограмме матрично-векторного умножения делалась копия входного вектора.
Таблица 5
Время решения задачи параллельной версией метода COCR(MRS)
V , Гц Число ядер
1 2 3 4
0.000584 13 мин 14.02 с 8 мин 22.28 с 6 мин 37.14 с 6 мин 14.12 с
0.001166 10 мин 43.45 с 6 мин 40.27 с 5 мин 19.34 с 5 мин 4.06 с
0.002326 8 мин 33.49 с 5 мин 23.37 с 4 мин 14.36 с 4 мин 5.35 с
0.004642 6 мин 41.05 с 4 мин 11.49 с 3 мин 18.06 с 3 мин 11.30 с
0.009261 5 мин 20.42 с 3 мин 22.34 с 2 мин 38.42 с 2 мин 37.18 с
0.01848 4 мин 9.73 с 2 мин 36.45 с 2 мин 4.53 с 2 мин 0.68 с
0.03687 3 мин 20.92 с 2 мин 6.34 с 1 мин 39.94 с 1 мин 35.81 с
0.07356 2 мин 41.95 с 1 мин 41.01 с 1 мин 20.88 с 1 мин 18.77 с
0.1468 2 мин 6.08 с 1 мин 19.57 с 1 мин 4.63 с 1 мин 1.83 с
0.2929 1 мин 40.92 с 1 мин 3.27 с 50.61 с 49.15 с
0.5843 1 мин 11.2 с 44.53 с 35.82 с 34.57 с
1.166 48.98 с 31.12 с 24.44 с 23.82 с
2.326 36.08 с 23.27 с 18.64 с 17.49 с
4.642 24.18 с 15.38 с 11.83 с 11.71 с
9.261 18.14 с 11.25 с 8.97 с 8.42 с
36.87 8.39 с 5.55 с 4.24 с 3.92 с
146.8 3.69 с 2.11 с 1.69 с 1.56 с
316.2 2.74 с 1.51 с 1.22 с 1.14 с
ЗАКЛЮЧЕНИЕ
Были рассмотрены вопросы, связанные с распараллеливанием в общей памяти итерационных методов решения конечно-элементных систем уравнений, получаемых при моделировании трехмерных магнитотеллурических полей с использованием векторного МКЭ. Вычислительные эксперименты показали преимущество метода COCR в сравнении с методами GMRES, BCG и COCG при решении данной задачи, поэтому для распараллеливания был выбран метод COCR.
Была проанализирована эффективность распараллеливания матрично-векторных умножений в итерационных методах решения конечно-элементных СЛАУ при использовании различных форматов хранения разреженных конечно-элементных матриц. Вычислительные эксперименты показали, что при относительно небольшом числе процессоров (2, 3, 4) следует отдать предпочтение хранению конечно-элементных матриц в разреженном строчно-столбцовом формате.
СПИСОК ЛИТЕРАТУРЫ
[1] Домников П.А. Моделирование трехмерных магнитотеллурических полей векторным методом конечных элементов // Актуальные проблемы электронного приборостроения. АПЭП-2008.: мат-лы IX Междунар. конф. / Новосибирск, 24-26 сентября, 2008. - Новосибирск: Изд-во НГТУ, 2008. - Т. 6. - С. 93-98.
[2] Ильин В.П. Методы и технологии конечных элементов. - Новосибирск: Изд-во ИВМиМГ СО РАН, 2007.
[3] Соловейчик Ю.Г., Рояк М.Э., Персова М.Г. Метод конечных элементов для решения скалярных и векторных задач: учеб. пособие. - Новосибирск: Изд-во НГТУ, 2007. - (Сер. «Учебники НГТУ»).
[4] Clemens M., Weiland T. Comparison of Krylov-Type Methods for Complex Linear Systems Applied to HighVoltage Problems // IEEE Transact. on Magnetics. - 1998. - Vol. 34, No. 5. - P. 3335-3338.
[5] Saad Y. Iterative Methods for Sparse Linear Systems. - 2nd ed., SIAM. - Philadelphia, 2003.
[6] Sogabe T., Zhang S.-L. A COCR method for solving complex symmetric linear systems // J. of Comp. and Appl. Mathem. - 2007. - № 199. - Р. 297-303.
[7] Van der Vorst H.A., Melissen J.B.M. A Petrov-Galerkin type method for solving Ax=b, where A is symmetric complex // IEEE Transact. on Magnetics. - 1990. - Vol. 26, No. 2. - P. 706-708.
[8] Zhou L., Walker H. Residual smoothing techniques for iterative methods // SIAM J. Sci. Comp. - 1994. - Vol. 15, No. 2. - Р. 297-312.
Домников Петр Александрович, аспирант кафедры прикладной математики Новосибирского государственного технического университета. Основное направление научных исследований - конечно-элементное моделирование трехмерных геоэлектромагнитных полей в изотропных и анизотропных средах. Имеет 10 публикаций.
E-mail: p_domnikov@mail.ru
Персова Марина Геннадьевна, кандидат технических наук, доцент, доцент кафедры прикладной математики Новосибирского государственного технического университета. Основное направление научных исследований - конечно-элементное моделирование электромагнитных полей в задачах геоэлектрики и электромеханики. Имеет 78 публикаций, в том числе 1 монографию.
E-mail: persova@fpm.ami.nstu.ru
Соловейчик Юрий Григорьевич, доктор технических наук, профессор, заведующий кафедрой прикладной математики Новосибирского государственного технического университета. Основное направление научных исследований - конечно-элементное моделирование электромагнитных и тепловых полей. Имеет более 100 публикаций, в том числе 1 монографию.
E-mail: kpmt@fpm.ami.nstu.ru
Вагин Денис Владимирович, аспирант кафедры прикладной математики Новосибирского государственного технического университета. Основное направление научных исследований - разработка методов моделирования трехмерных геоэлектромагнитных полей на нерегулярных сетках. Имеет 5 публикаций.
E-mail: vdvwk@mail.ru
P.A. Domnikov, M.G. Persova, Yu.G. Soloveichik, D.V. Vagin
Vector finite element modeling of 3D magnetotellurics and its performance capabilities for parallelization on shared memory processors.
This article is devoted to the problem of 3D magnetotellurics modeling by vector finite element method. Variational formulation for vector finite element method with separation of 1D part of field is exploited. Comparison of iterative solvers such as BCG, COCG, COCR with minimal residual smoothing and GMRES for solving finite element systems, arising in 3D magnetotellurics problems is carried out. Parallelization of iterative solvers on shared memory processors is performed. Problems of load balancing for different sparse storage schemes of finite element matrices are considered.
Key words: finite element method, parallel programming, electromagnetic field.
