close

Вход

Забыли?

вход по аккаунту

?

Использование CUDA Blas в решении квантовой задачи трёх частиц методом конечных элементов высоких порядков.

код для вставкиСкачать
УДК 530.145.61, 519.632.4
Вестник СПбГУ. Сер. 4, 2010, вып. 2
Д. А. Пузырев, Е. А. Яревский
ИСПОЛЬЗОВАНИЕ CUDA BLAS В РЕШЕНИИ КВАНТОВОЙ ЗАДАЧИ
ТРЁХ ЧАСТИЦ МЕТОДОМ КОНЕЧНЫХ ЭЛЕМЕНТОВ
ВЫСОКИХ ПОРЯДКОВ∗
Введение. Использование графических процессоров (GPU) для высокопроизводительных вычислений привлекает в последнее время все большее внимание научного
сообщества. Специфика архитектуры графических процессоров делает их особенно эффективными в вычислительных задачах, параллельных по данным. Одной из реализаций технологии вычислений общего назначения на GPU (GPGPU) является NVIDIA
CUDA. Использование технологии CUDA в некоторых приложениях может позволить
ускорить вычисления на порядки по сравнению с центральным процессором (CPU)
и добиться производительности, сравнимой с производительностью суперкомпьютеров.
Одной из областей вычислительной физики, требующей использования параллельных вычислений или средств аппаратного ускорения, является численное решение малочастичных квантовых задач. В данной работе исследовалась эффективность использования библиотеки CUDA BLAS для решения квантовой задачи трёх частиц, описываемой уравнением Шрёдингера.
Теоретическое описание проблемы. Полная волновая функция квантовой трёхчастичной системы с нулевым угловым моментом может быть описана уравнением
Шрёдингера в трёхмерных координатах. В качестве таких координат в данной работе
выбраны координаты Якоби: x – расстояние между двумя частицами, y – расстояние
между центром масс пары и третьей частицей, а φ – угол между векторами x и y. Уравнение Шрёдингера для дискретного спектра системы z в этих координатах имеет вид:
∂2
1 ∂2
y
−
x−
−
μ1,23 y ∂y 2
μ23 x ∂x2
1
1
∂2
∂
+ cot φ
+
∂φ2
∂φ
+ V (x, y, φ) ψ(x, y, φ) = zψ(x, y, φ). (1)
1
+
μ1,23 y 2
μ23 x2
Приведённые массы μ определены в терминах масс частиц m. Полный потенциал
V (x, y, φ) представляет собой сумму трёх парных взаимодействий. Волновая функция
ψ(x, y, φ) – ограниченная функция своих аргументов.
Численный метод. Уравнение (1) решалось с помощью метода конечных элементов [1, 2]. Применение метода конечных элементов (МКЭ) к задачам малочастичной
квантовой физики заметно отличается от большинства инженерных приложений этого
метода. Дело в том, что требуемая относительная точность решения задачи довольно
высока и должна быть не хуже чем 10−5 –10−6 . Достичь такой точности существенно
проще, применяя элементы высоких порядков.
В нашей реализации МКЭ трёхмерное пространство в координатах x, y и c = cos φ
разбивается на некоторое количество областей, имеющих форму параллелепипеда
и пронумерованных индексом i. На каждом элементе выбирается некоторое количество элементных функций, пронумерованных индексом m. Полная волновая функция
∗ По материалам доклада на юбилейном семинаре «Вычислительная физика» 29–30 октября
2009 г., С.-Петербург.
© Д. А. Пузырев, Е. А. Яревский, 2010
123
представляется следующим образом:
ψ(x, y, c) =
i
vim fim (x, y, c).
m
Функции fim (x, y, c) в методе конечных элементов обладают следующим свойством:
fim (x, y, c) ≡ 0 для (x, y, c) ∈
/ элемент i.
Мы выбираем элементные функции fim (x, y, c) в виде произведения одномерных функций по каждой координате. Вариационные коэффициенты vim и энергия z получаются
в результате решения вариационного уравнения. Приведённое выше разложение волновой функции приводит к следующей обобщённой задаче на собственные значения:
H̃v = z S̃v ,
(2)
где
(H̃)im,jk = fim |H|fjk , (S̃)im,jk = fim |fjk ,
а гамильтониан H определён в уравнении (1).
Вычислительные проблемы. Численное решение задачи (2) распадается на три
этапа: вычисление матричных элементов (H̃)im,jk , решение линейных систем с относительно разреженными матрицами и решение обобщённой задачи на собственные значения для матриц больших размерностей.
Соотношение требуемых вычислительных ресурсов для каждого из этих этапов зависит от требуемой точности расчётов и от изучаемой физической системы. Для решаемых в рамках развиваемого подхода задач (см., например, работы [2, 3], затраты
ресурсов на последний этап достаточно малы, и основное время тратится на расчёт
матричных элементов (H̃)im,jk и решение линейных систем.
В настоящей работе мы фокусируем наши усилия на ускорении вычисления матричных элементов. В силу специального, факторизованного, выбора элементных функций,
вычисление матричных элементов кинетической части гамильтониана (1) не представляет проблем и сводится к вычислению трёх одномерных интегралов. Однако вычисление матричных элементов потенциала сводится к трёхмерному интегралу
(3)
dx dy dc fim (x, y, c)V (x, y, c)fjk (x, y, c),
который не может быть представлен как произведение одномерных интегралов. Интегрирование по координатам x и y производится по квадратурным формулам Гаусса
с числом узлов, равным нескольким десяткам. Количество таких интегралов оценива2
ется как ∼ NelemNfunc
/2 и может быть весьма велико, так как количество элементных
функций Nfunc для применяемого нами МКЭ высоких порядков достигает нескольких
сотен. Таким образом, именно вычисление интегралов (3) составляет основную вычислительную задачу, которая будет решаться с использованием средств GPU.
Вычисление интегралов сводится к последовательности перемножения матриц относительно небольших размерностей (порядка нескольких десятков). В NVIDIA CUDA
включена библиотека CUBLAS [4], в которой реализованы все основные функции Basic
Linear Algebra Subprograms (BLAS), включая операцию перемножения матриц. Возможно использование CUBLAS в двух режимах: «thunking» и «non-thunking». В режиме «thunking» происходит обращение к функциям CUBLAS без изменения исходного
124
10- 2
CUBLAS Non-thunking
CUBLAS Thunking
MKL
Время
10- 3
10- 4
10- 5
0
50
100
150
M
200
250
300
Рис. 1. Время перемножения матриц размера M × M для двух режимов CUBLAS и Intel MKL BLAS
кода программы, при этом при каждом обращении к функции нужные данные из оперативной памяти автоматически копируются в память графического процессора, после
чего производятся вычисления и результаты расчёта копируются обратно. В режиме
«non-thunking» управление выделением памяти и копированием данных на GPU и обратно в оперативную память производится с помощью специальных служебных функций, а вызов самой функции CUBLAS запускает расчёт на графическом процессоре.
Известно, что в случае матриц больших размерностей (порядка нескольких тысяч)
использование CUBLAS в любом режиме приводит к уменьшению времени расчёта
на один-два порядка. Перемножение матриц небольшого размера при использовании
«thunking»-вызовов CUBLAS неэффективно по сравнению с CPU (рис. 1), так как
в этом случае основное время работы процедуры тратится на относительно медленные операции выделения памяти и копирования данных.
При использовании «non-thunking»-вызовов, как видно на рисунке, даже на небольших матрицах можно добиться ускорения расчёта по сравнению с CPU. В данном случае на каждые 100 операций умножения приходилось две операции копирования данных в память графического процессора и обратно в оперативную память, что примерно
соответствует соотношению количества операций в реальных расчётах для решаемой
задачи. При расчёте на CPU использовалась библиотека Intel MKL BLAS.
Результаты и выводы. Для исследования эффективности работы GPU при расчётах реальных квантовомеханических систем нами были проведены расчёты тримера
гелия 4 He3 [3, 5]. Использовались «non-thunking»-вызовы CUBLAS, при этом число операций выделения памяти графического процессора и копирования данных было сведено
к необходимому минимуму. Время расчёта матричных элементов для различного числа
узлов интегрирования приведено на рис. 2. Как видно, на данном этапе скорость расчёта матричных элементов на GPU приближена к CPU, при увеличении числа узлов интегрирования эффективность GPU начинает расти. Дальнейшее развитие оборудования
125
Время
500
CUBLAS
MKL
400
300
200
100
0
10
20
30
40
50
60
x,y
Nelem
70
80
Рис. 2. Время расчёта матричных элементов на
GPU и CPU
и технологий GPGPU, а также работа над переносом программы на GPU позволяет
ожидать значительного опережения скорости работы по сравнению с CPU.
Литература
1. Zienkiewicz O. C., Taylor R. L. The finite element method. The basis. Vol. 1. ButterworthHeinemann, 2000.
2. Elander N., Yarevsky E. Finite-element calculations of the antiprotonic helium atom including
relativistic and QED corrections // Phys. Rev. (A). 1997. Vol. 56. N 3. P. 1855–1864.
3. Salci M., Yarevsky E., Levin S. B., Elander N. Finite element investigation of the ground
states of the helium trimers 4 He3 and 4 He2 –3 He // Int. J. Quant. Chem. 2007. Vol. 107. P. 464–468.
4. CUDA BLAS library documentation, http://www.nvidia.com/cuda.
5. Motovilov A. K., Sandhas W., Sofianos S. A., Kolganova E. A. Binding energies and scattering observables in the 4 He3 atomic system // European Phys. J. (D). 2001. Vol. 13. P. 33–41.
Статья поступила в редакцию 22 декабря 2009 г.
Документ
Категория
Без категории
Просмотров
4
Размер файла
263 Кб
Теги
методов, blast, квантовое, использование, элементов, трех, высокие, задачи, решение, конечный, cuda, порядков, части
1/--страниц
Пожаловаться на содержимое документа