Метод фурье для уравнения пуассона

Видео:7.2 Задача 1. Краевая задача для уравнения ПуассонаСкачать

Аннотация научной статьи по математике, автор научной работы — Мингалев Олег Викторович, Мельник Михаил Николаевич

Для численного решения краевых задач для уравнения Пуассона в 2-мерном прямоугольнике и 3-мерном параллепипеде на регулярной сетке с большим числом узлов методом быстрого преобразования Фурье (БФП) разработаны несколько новых приемов, которые позволяют эффективно использовать параллельные вычисления как на ядрах центрального процессора, так и на графических процессорах (GPU). Создан набор программ для случая периодических граничных условий, а также для случаев однородных граничных условий Дирихле и Неймана. Эти программы дают решение с 4-м порядком точности и свободны от ограничения на число шагов сетки по каждому измерению в исходном методе БФП. Программы имеют простой и удобный интерфейс, а также максимально возможный уровень параллельности, и позволяют достаточно быстро решать указанные задачи с числом узлов сетки порядка 109 и более.

Видео:УМФ, 01.12, решение задач Лапласа и Пуассона в случае неоднородных граничных условийСкачать

Похожие темы научных работ по математике , автор научной работы — Мингалев Олег Викторович, Мельник Михаил Николаевич

Видео:9. Уравнение ПуассонаСкачать

NUMERICAL SOLUTION OF BOUNDARY PROBLEMS FOR THE EQUATIONS OF POISSON BY FAST FOURIER TRANSFORM USING PARALLEL CALCULATIONS

For the numerical solution of boundary value problems for the Poisson equation in a 2D and a 3D on a regular grid with a large number of nodes, the Fast Fourier Transform (FFT) method has developed several new techniques that make it possible to efficiently use parallel computing on both the cores CPU and on graphics processing units (GPU). A set of programs has been created for the case of periodic boundary conditions, as well as for the cases of homogeneous Dirichlet and Neumann boundary conditions. These programs have a solution 4th order of accuracy and haven’t the restriction on the number of grid steps for each dimension in the original FFT method. Programs have a simple and usable interface, as well as the highest possible level of parallelism, and allow you to quickly solve these problems with the number of grid nodes of the order of 109 or more.

Видео:7.6 Задача 5. Краевая задача для уравнения ПуассонаСкачать

Текст научной работы на тему «Численное решение краевых задач для уравнения Пуассона методом быстрого преобразования Фурье с использованием параллельных вычислений»

DOI: 10.25702/^^2307-5252.2018.9.5.165-182 УДК 533.95 + 519.63

О. В. Мингалев, М. Н. Мельник

ЧИСЛЕННОЕ РЕШЕНИЕ КРАЕВЫХ ЗАДАЧ ДЛЯ УРАВНЕНИЯ ПУАССОНА МЕТОДОМ БЫСТРОГО ПРЕОБРАЗОВАНИЯ ФУРЬЕ С ИСПОЛЬЗОВАНИЕМ ПАРАЛЛЕЛЬНЫХ ВЫЧИСЛЕНИЙ

Для численного решения краевых задач для уравнения Пуассона в 2-мерном прямоугольнике и 3-мерном параллепипеде на регулярной сетке с большим числом узлов методом быстрого преобразования Фурье (БФП) разработаны несколько новых приемов, которые позволяют эффективно использовать параллельные вычисления как на ядрах центрального процессора, так и на графических процессорах (GPU). Создан набор программ для случая периодических граничных условий, а также для случаев однородных граничных условий Дирихле и Неймана. Эти программы дают решение с 4-м порядком точности и свободны от ограничения на число шагов сетки по каждому измерению в исходном методе БФП. Программы имеют простой и удобный интерфейс, а также максимально возможный уровень параллельности, и позволяют достаточно быстро решать указанные задачи с числом узлов сетки порядка 109 и более.

уравнения Пуассона, краевая задача, метод быстрого преобразования Фурье, параллельные вычисления.

O. V. Mingalev, M. N. Melnik

NUMERICAL SOLUTION OF BOUNDARY PROBLEMS FOR THE EQUATIONS OF POISSON BY FAST FOURIER TRANSFORM USING PARALLEL CALCULATIONS

For the numerical solution of boundary value problems for the Poisson equation in a 2D and a 3D on a regular grid with a large number of nodes, the Fast Fourier Transform (FFT) method has developed several new techniques that make it possible to efficiently use parallel computing on both the cores CPU and on graphics processing units (GPU). A set of programs has been created for the case of periodic boundary conditions, as well as for the cases of homogeneous Dirichlet and Neumann boundary conditions. These programs have a solution 4th order of accuracy and haven’t the restriction on the number of grid steps for each dimension in the original FFT method. Programs have a simple and usable interface, as well as the highest possible level of parallelism, and allow you to quickly solve these problems with the number of grid nodes of the order of 109 or more.

Poisson’s equation, of boundary value problems, the Fast Fourier Transform, Parallel calculation.

Проблема быстрого численного решения краевых задач для уравнения Пуассона (см., например, [1]) в 2-мерном прямоугольнике и 3-мерном

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

Для этих целей возможно использовать только такие численные методы, в которых основной объем вычислений может выполняться в параллельном режиме на графических процессорах (GPU). К числу таких относится метод быстрого преобразования Фурье (далее БФП) (см., например, 2), который позволяет решать задачу с периодическими граничными условиями, а также задачи с однородными граничными условиями Дирихле и Неймана.

Отметим, что в платные библиотеки MKL и IMSL языка FORTRAN входят подпрограммы решения краевых задач для уравнения Пуассона в 2-мерном прямоугольнике и 3-мерном параллепипеде методом БФП, однако они имеют ряд недостатков, которые делают их непригодными для достижения указанных выше целей.

Метод быстрого преобразования Фурье допускает два варианта: со 2-м и 4-м порядком точности, и основан на алгоритме Кули-Туки (Кули-Тьюки) дискретного быстрого преобразования Фурье периодического 1 -мерного массива

с размерностью N = 2m , что дает соответствующее ограничение на числа Nxq ,

Nyо, Nzq шагов сетки по каждому измерению X, Y, Z вида

Nxo = 2mx, Nyo = 2my, Nzo = 2mz , (1)

где mx, my, mz e N — натуральные числа. Из этого ограничения вытекает условие

на соотношения размеров параллепипеда (прямоугольника) Lx , Ly , Lz по соответствующим измерениям X, Y, Z вида

Ly/Lx = 2My,х , Lz/Lx = 2Mz,x , My, e Z, M, e Z , (2)

где My, = my — mx, Mz, = mz — mx e Z — целые числа.

Это ограничение весьма неудобно для построения численных моделей крупномасштабных процессов в космической плазме. Отметим, что алгоритм Кули-Туки может быть обобщен также и для случая, когда размерность массива N имеет вид произведения степеней простых чисел 2, 3, 5, 7, то есть

дг _ 2>»2.3»% . 5Щ . i™! ^ где т2,т3,т5,т1 е Z+ — целые неотрицательные

числа. Однако в этом случае он существенно усложняется.

Для численного решения задачи Дирихле и задачи Неймана можно использовать два способа. Первый способ состоит в использовании соответственно дискретного

синус или косинус-преобразования Фурье. Второй способ состоит в сведении исходной задачи к периодической задаче в прямоугольной области с удвоенными размерами путем продолжения правой части специальным образом. Отметим, что быстрое дискретное синус и косинус-преобразование Фурье периодического 1 -мерного массива с размерностью N сводятся к обычному быстрому дискретному преобразованию Фурье периодического 1 -мерного массива с удвоенной размерностью 2N . При этом программы синус и косинус-преобразования существенно сложнее исходного преобразования Фурье. Кроме того, при таком сведении задачи Неймана к периодической автоматически достигается аппроксимация однородного граничного условия Неймана с тем же порядком точности, что и аппроксимация уравнения Пуассона. Поэтому второй способ для наших целей является более удобным, и в работе предложен прием сведения задачи Дирихле или Неймана с однородными граничными условиями в прямоугольнике [0;Ьх] х [0;Ьу ] или параллепипеде

[0; Ьх ] х [0; Ьу ] х [0; Ь2 ] к периодической задаче в аналогичной области с

удвоенными размерами по каждому измерению. При этом автоматически получается конечноразностная аппроксимация граничного условия Неймана того же порядка точности, что и аппроксимация уравнения Пуассона, в частности, аппроксимация 4-го порядка точности. Это является важным новым методическим результатом.

Также в работе предложен новый прием, который позволяет снять ограничение (1) на число шагов сетки по каждому измерению, то есть оставляет только условие (2), и состоит в следующем. Если условие (1) нарушено, то по исходной сетке с шагом й0 создается вспомогательная сетка с более мелким шагом И и числами шагов сетки по каждому измерению Nx, N у , N2 , которые удовлетворяют условию (1). При этом

числа тх, Шу, тг определяются через исходные числа шагов Nxo , N у 0 , N20 по формулам

Шх = N2 Nx0 ] +1> Шу =[^2 Nyo ] +1, т2 = [1о§2 Nzo ] +1, (3)

где через [а] обозначена целая часть действительного числа а. Для получения

дискретной задачи на вспомогательной сетке выполняем интерполяцию правой части уравнения с исходной сетки на вспомогательную сетку при помощи интерполяционного полинома Лагранжа с порядком на 2 большим используемого порядка аппроксимации уравнения. То есть в случае метода 2-го порядка точности используем полином Лагранжа 4-го порядка — аппроксимация по 4-м ближайшим узлам на измерение, а в случае метода 4-го порядка точности используем полином Лагранжа 6-го порядка — аппроксимация по 6-ти ближайшим узлам на измерение. После этого методом БФП решается задача на вспомогательной сетке, а потом полученное численное решение интерполируется обратно на исходную сетку полностью аналогично описанной интерполяции. Многочисленные тестовые расчеты показали высокую точность и эффективность этого приема.

В работе рассматриваются основные детали программы, которая выполняет следующие шаги.

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

2) Построение удовлетворяющей условию (1) вспомогательной сетки с меньшим шагом И , и интерполяция правой части уравнения на эту сетку, которую удобно выполнять в параллельном режиме.

3) Вычисление собственных чисел 1-мерной дискретной периодической задачи по каждому измерению на вспомогательной сетке.

4) Вычисление массива коэффициентов дискретного преобразования Фурье правой части на вспомогательной сетке. Эта первая из двух главных частей программы, которая содержит почти половину всех вычислений. Она выполняется в результате последовательного выполнения по каждому измерению 1-мерного быстрого преобразования Фурье для всех 1-мерных сечений массива правой части вдоль данного измерения. В ходе последнего процесса каждое 1 -мерное сечение массива правой части вдоль данного измерения независимо обрабатывается подпрограммой 1-мерного БФП. Поэтому эта часть программы удобна для выполнения в параллельном режиме как на ядрах центрального процессора, так и на графическом процессоре.

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

6) Вычисление сеточного массива решения при помощи обратного дискретного преобразования Фурье. Эта вторая из двух главных частей программы, которая по объему вычислений и методике выполнения аналогична прямому преобразованию Фурье, описанному в пункте 4).

7) Обратная интерполяция решения со вспомогательной сетки на исходную сетку, аналогичная прямой интерполяции в пункте 2).

8) Численное дифференцирование решения для нахождения соответствующего векторного поля, использующегося в модели.

Отметим, что все описанные выше шаги удобны для выполнения в параллельном режиме.

1. Сведение задач Дирихле и Неймана к периодической задаче.

Обозначим в пространствах К2 и М3 векторы ортонормированного базиса декартовой системы координат через ех , еу и е2, а векторы в этих

пространствах будем рассматривать в виде:

2 ^ л- =хех + уеу еК , х=хех+ уеу+ .

Будем рассматривать 2-мерный открытый исходный прямоугольник П21 и

прямоугольник П2 2 с удвоенными размерами, которые определяются

а также аналогичные 3-мерные параллепипеды исходный П3 ^ и «удвоенный» П3 2 , которые определяются формулами

Пзд =(0;Ьх)х(0;Ьу)х(0;Ь2), =(0;2!х)х(0;2Ьу)х(0;2Ьг) . (5)

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

Рассмотрим краевые задачи Дирихле и Неймана с однородными граничными условиями и финитной правой частью для уравнения Пуассона

А и(х) = /(х) , х е Пи>1, (6)

где правая часть /(х) — заданная функция. Будем считать, что правая часть

У(х) е С^ПиД^ , то есть является достаточно гладкой (1 раз непрерывно

дифференцируема) на замкнутом прямоугольнике (параллепипеде) Пп,1 и финитной (равной нулю в некоторой окрестности границы ЗП^ ), где п = 2,3 -размерность пространства. Однородное граничное условие Дирихле имеет вид

и( х ) = 0, х е ЗПп1. (7)

Однородное граничное условие Неймана имеет вид

где н (х) — единичная внешняя нормаль к границе ЗП^ . Рассмотрим сведение задачи Дирихле (6), (7) или задачи Неймана (6), (8) с однородными граничными условиями в прямоугольнике или параллепипеде Пп 1 к соответствующей периодической краевой задаче в удвоенном прямоугольнике или параллепипеде Пп2 . Для этого необходимо соответствующим образом продолжить правую

часть /(х) с множества Пп,1 на удвоенный прямоугольник или параллепипед Пп,2 . Обозначим такое продолжение в случае задачи Дирихле через FD(х), а в случае задачи Неймана через FN( х) . Для продолжения FD(х) по каждой переменной используется нечетное продолжение относительно правого края Ьх , Ьу , Ьг исходного интервала. В 2-мерном случае продолжение FD(х) определяется формулами

^ ( х у ) = /( х, у) , ^ ( х + Ьх, уУ) = -/(Ьх- х, у ) , ^ ( ^ у + Ьу ) = -/( х, Ьу- у) , ^ ( х + Ьх, у + Ьу ) = /(Ьх- х, Ьу- у ),

( х, у ) е П2,1 = [0; Ьх ] х[0; Ьу ].

В 3-мерном случае продолжение FD( х) получается в результате последовательного применения следующих формул:

(х У, *) = у(х, У, *) , (х + Ьх, У, * ) = -/(Ьх- х, У, * ) Рв (x, У + Ьу, * ) = -/(х, Ьу- У, * ) ,

Рв (х + Ьх, У + ЬУ, *) =1 (Ьх — х, ЬУ — У, *)

при (х, У, * )е Пз,1 = [0; Ьх ] х [0; Ьу ] х [0; Ь* ], Рв ( x, У, * + Ь* ) = -/ ( х, У, Ь*- * )

при (х,У,*) е [0;2Ьх] х[0;2Ьу] х[0;Ь*].

Для продолжения Ру(х) по каждой переменной используется четное

продолжение относительно правого края исходного интервала. В 2-мерном случае оно о пределяется формулами

Рм (х у) = /(х, У) , Рм (х + Ьх, У) = у(Ьх- х, У) ,

Рм (х, У + Ьу ) = /(х, Ьу- у) , (х + Ьх, у + Ьу ) = /(Ьх- х, Ьу- у) , I (11)

(х, у) е П2,1 = [0; Ьх ] х [0; Ьу ], ^

а в 3-мерном случае получается в результате последовательного применения аналогичных формул:

Рм (х У, *) = У (х, У, *), Рм (х + Ьх, У, *) = У (Ьх- х, У, *), Рм ( x, У + Ьу, * ) = У ( х, Ьу- У, * ),

Рм (х + Ьх, У + Ьу, *) =1 (Ьх- х, Ьу- У, * ) при (х, у, * )е Пз,1 = [0; Ьх ] х [0; Ьу ] х [0; ^ ],

Рм ( х У, * + Ь* ) = у ( х, У, Ь*- * ) при ( х, у, * )е [0;2Ьх ] х [0;2Ьу ] х [0; Ь* ].

Рассмотрим в удвоенном прямоугольнике или параллепипеде Пи 2 периодические краевые задачи для уравнения Пуассона

ди (х ) = Р (х), х е Пи,2, (13)

с продолженными правыми частями Р( х) = х) и Р( х) = Рм( х) . Используя

представление решения его рядом Фурье по собственным функциям, можно доказать следующее утверждение.

Утверждение. Пусть /(х)еС^Пид) , где область Пи1 определяется формулами (4) и (5), а функции х) и Рм(х) являются продолжениями функции /(х) по формулам (9), (11) и (10), (12) соответственно на

определяемую формулами (4) и (5) удвоенную область Пи2 . Если ив (х) является классическим решением периодической краевой задачи (13) при

F (x) = Fd(x) , то его ограничение на Пид является классическим решением задачи Дирихле (6), (7). Если UN(x) является классическим решением периодической краевой задачи (13) при F (x) = FN(x) , то его ограничение на Пи,1 является классическим решением задачи Неймана (6), (8).

2. Интерполяция на вспомогательную сетку и обратно.

Рассмотрим вопрос построения вспомогательной сетки. Будем считать, что

размеры параллепипеда (прямоугольника) Lx , Ly , Lz удовлетворяют условию (2), и в задаче из физических условий введена сетка с шагом h = Lx/Nxq = LyjNyft = Lz/Nzq . Если числа шагов сетки по каждому

измерению Nx0 , Ny0 , Nz0 не удовлетворяют условию (1), то вводим

вспомогательную сетку с меньшим шагом h = Lx¡Nx = LyjNy = Lz/Nz , для

которой новые числа шагов сетки по каждому измерению Nx , N y , Nz

удовлетворяют условию (1). При этом числа mx, m.y, mz определяются через

исходные числа шагов Nx0 , Ny0, Nz0 по формулам (3). При помощи целых

векторов (мультииндексов) к = кх ex + ky ey и к = кх ex + ку ey при n = 2 , и

к = кх ex + ky ey + kz ez и к = kx ex + ky ey + kz ez при n = 3 узлы этих сеток обозначим как

r(k) = ho-k, f(k) = hk, к,к&Ъп, п = 2,3, (14)

Будем считать, что в исходной физической задаче задан массив значений Fh (к) = F(r (к)) правой части F(x) в узлах исходной сетки на замыкании

«удвоенной» области ПИ;2 . Для получения массива значений Fh (к^ = F^r (к^

правой части в узлах вспомогательной сетки выполняем интерполяцию с исходной сетки на вспомогательную при помощи интерполяционного полинома Лагранжа 6-го порядка — аппроксимация по 6-ти ближайшим узлам на измерение, то есть по 36 узлам при n = 2 и по 216 узлам при n = 3 . Применительно к рассматриваемому случаю эту формулу интерполяции при n = 2 удобно представить в форме

= Z Z Fhilhklh0] + Py

U S’x, Px + (px — ‘x ) X

где через [Я ] и <Я>= Я — [Я] обозначены целая и дробная части действительного числа Я, через ■ к¡И<)^ = - кх /+ - ку /

р = рхех + ру еу е Z2 обозначены соответствующие векторы с целыми координатами, и через

обозначен символ Кронекера. При п = 3, используя аналогичные обозначения, формулу интерполяции удобно представить в форме

🌟 Видео

7.4 Задача 3. Краевая задача для уравнения ПуассонаСкачать

29. Адиабатический процесс. Уравнение ПуассонаСкачать

7.3 Задача 2. Краевая задача для уравнения ПуассонаСкачать

7.5 Задача 4. Краевая задача для уравнения ПуассонаСкачать

7.9 Задача 8. Краевая задача для уравнения ПуассонаСкачать

Уравнения математической физики 15+16 Задача Дирихле для уравнения Лапласа - Пуассона в кругеСкачать

Уравнения математической физики 11 Формула Пуассона для уравнения теплопроводностиСкачать

Уравнение колебаний струны. Метод разделения переменных. Метод ФурьеСкачать

Билет №04 "Потенциал электростатического поля"Скачать

Практическое занятие. Численное решение уравнений Лапласа и ПуассонаСкачать

Метод Фурье для неоднородного уравнения теплопроводностиСкачать

Уравнение ПуассонаСкачать

7.8 Задача 7. Краевая задача для уравнения ПуассонаСкачать

Задача Дирихле для круга. Уравнение ЛапласаСкачать