Параллельные алгоритмы для исследования системы длинных джозефсоновских переходов

Параллельные алгоритмы для исследования системы длинных джозефсоновских переходов*

М.В. Башашин1,3, Е.В. Земляная1,3, М.И. Зуев1, А.В. Нечаевский1, И. Р. Рахмонов2, Д.В. Подгайный1, О.И. Стрельцова1,3, Ю.М. Шукринов2,3
1Лаборатория информационных технологий ОИЯИ
2Лаборатория теоретической физики им. Н.Н. Боголюбова ОИЯИ
3Государственный университет «Дубна»

*Исследование выполнено за счет гранта Российского научного фонда (проект № 18-71-10095).

В работе представлены результаты исследования эффективности параллельных реализаций вычислительной схемы для расчета вольт-амперных характеристик системы длинных джозефсоновских переходов [1,2]. В работах [3,4] дано описание разработанного параллельного алгоритма и представлены результаты тестирования MPI-версии расчета ВАХ в зависимости от вычислительных параметров. В настоящей работе представлены разработанные параллельные реализации на основе другого инструментария: OpenMP-реализация для проведения расчетов на системах с общей памятью, CUDA-реализация для проведения расчетов на графических процессорах NVIDIA. Разработка, отладка и профилирование параллельных приложений проводились на учебно-тестовом полигоне гетерогенной вычислительной платформы HybriLIT, а расчеты проводились на суперкомпьютере ГОВОРУН [2].

Постановка задачи. Рассматривается обобщенная модель, учитывающая индуктивную и емкостную связь между длинными джозефсоновскими переходами (ДДП) [1]. Система связанных ДДП предполагается состоящей из  сверхпроводящих слоев c промежуточными диэлектрическими слоями длины L.

Фазовая динамика системы N ДДП с учетом емкостной и индуктивной связи между контактами описывается начально-краевой задачей для системы дифференциальных уравнений относительно разности фаз \phi _{l}(x,t) и напряжения V _{l}(x,t) на каждом  l-ом контакте (l=1,2,...N) . В безразмерном виде система уравнений имеет вид [1]:

    \[ \left\{\begin{matrix} \frac{\partial \phi }{\partial t}& = C\cdot V, 0<x<L,t>0\\ \frac{\partial V }{\partial t}& =L^{_{-1}}\frac{\partial^2\phi }{\partial x^2}-\beta V-sin(\phi )+I(t); \end{matrix}\right. \]

    \[ \phi =\begin{pmatrix} \phi_{1}\\ \phi_{2}\\ \cdot \cdot \cdot \\ \phi_{N}\\ \end{pmatrix}, V =\begin{pmatrix} V_{1}\\ V_{2}\\ \cdot \cdot \cdot \\ V_{N}\\ \end{pmatrix}, \]

где  L — матрица индуктивной связи, С — матрица емкостной связи:

    \[ L=\begin{pmatrix} 1& S& 0& \cdots& \square& 0& S \\ \vdots& \vdots& \vdots& \vdots& \vdots& \vdots& \vdots&\\ \cdots& 0& S& 1& S& 0& \cdots&\\ \vdots& \vdots& \vdots& \vdots& \vdots& \vdots& \vdots&\\ \square& \square& \square& \square& \square& \square& \square& \\ S& 0& \cdots& \cdots& 0& S& 1& \end{pmatrix}, C=\begin{pmatrix} D^{C}& S^{C}& 0& \cdots& \cdots& 0& 0& S^{C}&\\ \vdots& \vdots& \vdots& \vdots& \vdots& \vdots& \square& \vdots&\\ \cdots& 0& S^{C}& D^{C}& S^{C}& 0& \cdots& \cdots&\\ \vdots& \vdots& \vdots& \vdots& \vdots& \vdots& \vdots& \vdots&\\ \cdots& \cdots& \cdots& \cdots& 0& S^{C}& D^{C}& S^{C}&\\ S^{C}& 0& \cdots& 0& S^{C}& 0& S^{C}& D^{C}& \end{pmatrix}, \]

β — параметр диссипации, параметр индуктивной связи S принимает значение на интервале 0<||<0.5; D c— эффективная электрическая толщина ДП, нормированная на толщину диэлектрического слоя; sc — параметр емкостной связи, l(t) — внешний ток.

Система уравнений дополняется нулевыми начальными условиями и граничными условиями:

    \[ V_{l}(0,t)=V_{l}(L,t)=\frac{\partial \phi (0,t)}{\partial x}=\frac{\partial \phi (L,t)}{\partial x}=0, l=1,2,...,N. \]

Также рассматривалась задача, когда граничные условия в направлении x задаются внешним магнитным полем [2].

При вычислении вольт-амперных характеристик (ВАХ) зависимость тока от времени выбиралась в виде ступенек (схематичное представление дано на рисунке ниже), т.е. решается задача при постоянном токе (I=I_{j}), а найденные  функции \phi_{l},V_{l}(l=1,2,...N) берутся в качестве начальных условий для решения зaдачи для тока &I_{j+1}&.

Вычислительная схема и схема распараллеливания. Для численного решения начально-краевой задачи строилась равномерная сетка по пространственной переменной (количество узлов сетки NX ). В системе уравнений (1) производная второго порядка по координате x аппроксимируется с помощью трехточечных конечно разностных формул на дискретной сетке с равномерным шагом \Delta x. Полученная система дифференциальных уравнений относительно значений \phi_{l},V_{l}(l=1,2,...N) в узлах дискретной сетки по x решается методом Рунге–Кутта четвертого порядка.

Для вычисления ВАХ проводится усреднение V_{l}(x,t) по координате и по времени. Для этого на каждом шаге по времени проводится интегрирование напряжения по координате методом Симпсона и усреднение

    \[ \bar{V_{l}}=\frac{1}{L}\int_{0}^{L}V_{l}(x,t)dx, \]

далее напряжение усредняется по времени с использованием формулы

    \[ <\bar{V_{l}}>=\frac{1}{T_{max}-T_{min}}\int_{T_{min}}^{T_{max}}\bar{V_{l}}(t)dt \]

и затем суммируется по всем ДП. Для интегрирования по времени используется метод прямоугольников.

Схема распараллеливания. При численном решении начально-краевой задачи методом Рунге-Кутта 4-го порядка по временной переменной, на каждом временном слое коэффициенты Рунге-Кутта (K_{i}) могут быть найдены независимо (параллельно) для всех NX узлов пространственной сетки и для всех N джозефсоноских переходов. При этом, коэффициенты K_{i}(i=1,2,3,4) определяются один за одним (последовательно). Таким образом, эффективно параллелизация проводится по NX \cdot N точкам. Также при проведении усреднений при расчетах ВАХ, вычисление интегралов может быть проведено в параллельном режиме.

Параллельные реализации. Для ускорения расчетов ВАХ были разработаны параллельные реализации описанной выше вычислительной схемы. Ниже приведены результаты исследования эффективности параллельных реализаций, проведенные при значениях параметров: L=10;I_{min}=0;I_(max)=1.1;\beta=0.2;N=1; T_(max)=200; \Delta t=0.04, количество узлов по пространственной переменной:

NX = 20048.

OpenMP-реализация.
Расчеты проводились:

  • на вычислительных узлах с процессорами Intel Xeon Phi 7290 (KNL: 16GB, 1.50 GHz, 72 core, поддерживается 4 потока на ядро — всего 288 логических ядер), компилятор Intel (Intel Cluster Studio 18.0.1 20171018);
  • на двухпроцессорных на вычислительных узлах с процессорами Intel Xeon Gold 6154 (Skylake: 24.75M Cache, 3.00 GHz, 18 core, поддерживается 2 потока на ядро — всего 72 логических ядра на узел).

Ниже приведены графики зависимости ускорения вычислений, получаемое при использовании параллельного алгоритма:

    \[ S=\frac{T_{1}}{T_{n}}; \]

T_{1} время вычислений с задействованием одного ядра, T_{n} — время вычислений на n-логических ядрах, от количества потоков (число которых равно количеству логических ядер) и график зависимости эффективности использования параллельным алгоритмом вычислительных ядер:

    \[ E=\frac{T_{1}}({n\cdot T_{n})}\cdot100 \% . \]

характеризующей масштабируемость параллельного алгоритма.

На Рис.1 и Рис. 2 представлены зависимости ускорения и эффективности от количества потоков  при проведении вычислений на узлах с процессорами Intel Xeon Phi, а на Рис.3 и Рис.4 ускорение и эффективность при подключении возможности использования инструкции AVX-512. Отметим, что использование этой инструкции позволило сократить время расчетов в среднем в 1.8 раза.

Рисунок 1
Рис. 1. График зависимости ускорения параллельных вычислений от количества потоков
Рисунок 2
Рис. 2. График зависимости эффективности использование параллельным алгоритмом вычислительных ядер от количества потоков
Рисунок 3
Рис. 3. График зависимости ускорения параллельных вычислений от количества потоков. Компиляция программы с опцией — xMIC-AVX512
Рисунок 4
Рис. 4. График зависимости эффективности использование параллельным алгоритмом вычислительных ядер от количества потоков. Компиляция программы с опцией — xMIC-AVX512

На Рис.5 и Рис. 6 представлены зависимости ускорения и эффективности от количества потоков  при проведении вычислений на узлах с процессорами Intel Xeon Gold 6154.

Рисунок 5
Рис. 5. График зависимости ускорения параллельных вычислений от количества потоков
Рисунок 6
Рис. 6. График зависимости эффективности использование параллельным алгоритмом вычислительных ядер от количества потоков

 

CUDA-реализация. Для расчетов на графических ускорителях от NVIDIA была разработана CUDA-реализация параллельного алгоритма. При вычислениях интегралов использовался алгоритм параллельной редукции с использованием shared memory.  Времена расчетов на графических ускорителях NvidiaTesla K40 и NvidiaTesla K80 представлены на диаграмме ниже.

Рис. 7. Время вычислений на GPU

Сравнительный анализ параллельных реализаций. В таблице 1 приведены результаты минимальных времен расчетов, полученных на различных вычислительных архитектурах.

Таблица 1. Наилучшие (минимальные) времена расчетов, полученные на различных вычислительных архитектурах (в минутах) и соответствующие ускорения вычислений.
Таблица 2. Времена расчетов (в минутах) на одном OpenMP потоке процессора Intel Xeon Gold 6154(Skylake) и графических ускорителях Nvidia Tesla K40 и Nvidia Tesla K80.

Наилучшее (минимальное) время для представленных выше значений параметров расчета ВАХ ДДП составляет 1290 минут и было достигнуто при вычислениях на двухпроцессорном узле с процессорами Intel Xeon Gold 6154.

[1] И. Р. Рахмонов, Ю. М. Шукринов, П. Х. Атанасова, Е. В. Земляная, М. В. Башашин, Влияние индуктивной и емкостной связи на вольтамперную характеристику и электромагнитное излучение системы длинных джозефсоновских переходов, Журнал Экспериментальной и Теоретической Физики, Том 151, № 1, 151-159, (2017)

[2]  Gh. Adam, M. Bashashin, D. Belyakov, M. Kirakosyan, M. Matveev, D. Podgainy, T. Sapozhnikova, O. Streltsova, Sh. Torosyan, M. Vala, L. Valova, A. Vorontsov, T. Zaikina, E. Zemlyanaya, M. Zuev.  IT-ecosystem of the HybriLIT heterogeneous platform for high-performance computing and training of IT-specialists. CEUR Workshop Proceedings (CEUR-WS.org) (2018 г.)

[3] М. В. Башашин, Е. В. Земляная, И.Р. Рахмонов, Ю.М. Шукринов, П. Х. Атанасова, А. В. Волохова, Вычислительная схема и параллельная реализация для моделирования системы длинных джозефсоновских переходов, Компьютерные исследования и моделирование (Computer Research and Modeling), Том 8, № 4, 593–604, (2016)

[4] E. V. Zemlyanaya, M. V. Bashashin, I. R. Rahmonov, Yu. M. Shukrinov, P. Kh. Atanasova, and A. V. Volokhova, Model of stacked long Josephson junctions: Parallel algorithm and numerical results in case of weak coupling, AIP Conference Proceedings, vol 1773, 110018, (2016)

 

Скачать файл с текстом