1 - УрФУ (Екатеринбург) 2 - ИММ УрО РАН (Екатеринбург) 3 - Институт иммунологии и физиологии (Екатеринбург) Исследование функциональных характеристик миокарда является сложной задачей и требует математических моделей, описывающих процессы на разных уровнях (клетка-ткань-орган). Целью работы является определение возможности применения библиотек автоматизированных научных вычислений для моделирования распространения потенциала действия с помощью параллельных вычислительных систем на примере модели Алиева-Панфилова. Представлен обзор существующих программных комплексов (OpenFOAM, OpenCMISS, FEniCS, Chaste), а также результаты тестирования производительности и масштабируемости реализации на основе FEniCS.
начиналась с учетом особенностей одной технологии, которая через три года успевает устареть, и появляется необходимость в переписывании значительной части кодовой базы проекта, но уже с принятием во внимание специфики новой, более совершенной технологии.
В данных условиях постоянно приходится адаптировать и оптимизировать программный код под развивающиеся технологии. Дальнейшей целью нашей работы является многоуровневое моделирование, что влечет за собой дополнительные сложности в оптимизации. Изменения какого-либо участка кода на нижнем уровне требуют изменений на вышестоящих уровнях. Чем больше оптимизаций, тем сильнее программный код отдаляется от математической нотации. ¾Читабельность¿ кода теряется для математика (и тем более для физиологов), не знакомого с внутренним строением современных вычислительных архитектур и техниками оптимизации. Программа становится сложной и запутанной, а дальнейшее развитие и сопровождение проекта – нетривиальным.
В связи с вышесказанным можно выделить следующие основные сложности, которые возникают в крупных научных проектах, требующих массивных вычислений: 1) сохранение близости программного кода к математическим формулировкам; 2) достижение высокой скорости параллельных вычислений.
Как правило, чем код ближе к математической нотации, тем хуже он оптимизирован и наоборот. Существует ли решение данных проблем одновременно? Мы задались этим вопросом в начале нашего проекта и изучили ряд систем, которые предлагают это решение: OpenFOAM, OpenCMISS, Chaste, FEniCS.
В данной статье проведен обзор данных систем, выбрана наиболее подходящая для выполнения проекта, описана модель левого желудочка сердца, которая используется для тестовых расчетов, сформулирована задача моделирования и представлены результаты тестирования на вычислительном кластере. 2
Рассчитываемая модель
Модель Алиева–Панфилова принадлежит к группе так называемых феномологических моделей. От
других видов эту группу отделяет наличие притягательного свойства простоты реализации, но в то же
время они неспособны воспроизвести количественно важные характеристики исследуемой системы. Тем
не менее, феномологические модели возбудимых сред предоставляют хорошее качественное соответствие
свойств объекта моделирования и полезны, например, при исследовании сердечной аритмии [
Система уравнений Алиева–Панфилова [
1v 2 + u ku(u a)(u 1) uv; ( v ku (u a 1)) ; (1) (2) где u – безразмерный трансмембранный потенциал клетки, v – безразмерная проводимость для калиевого тока, dij – коэффициенты диффузии, а 0, a, k, 1, 2 – параметры модели клетки. Везде индексы i; j принимают значения 1; 2; 3. Первое слагаемое в правой части (1) отвечает за описание пространственного изменения потенциала. Другие слагаемые моделируют быстрые процессы изменения трансмембранного потенциала, в то время как уравнение (2) описывает медленные процессы.
Для целей тестирования вычислительных инструментов мы использовали симметричную модель формы
левого желудочка сердечной мышцы [
Модель Алиева–Панфилова была выбрана нами ввиду того, что, несмотря на свою относительную простоту, она является комбинацией нелинейного уравнения в частных производных параболического типа и жестких обыкновенных дифференциальных уравнений и отражает основные математические особенности моделей, которые используются для исследования свойств миокарда. Существующие программные средства моделирования свойств тел сложной формы 3.1
OpenFOAM может быть охарактеризована как платформа для исследования моделей турбулентных
течений, транспортных уравнений и других типов моделей, связанных с перераспределением различных
веществ в пространстве. В её основе лежит метод конечных объемов, который наиболее пригоден для
решения уравнений в частных производных, описывающих явления переноса. OpenFOAM состоит из двух
частей. Одна из них представляет собой набор классов, которые используются в реализации численных
процедур. Другая часть – библиотека программ-решателей ("solvers"). Такое разделение предоставляет
программную среду для моделирования с помощью высокоуровнего кода, который напрямую отражает
интегрируемые уравнения [
fvm : : l a p l a c i a n (D_T,T ) ) ; 3.2
OpenCMISS – это часть международного проекта Physiome [
Chaste 3.4 FEniCS
Другим пакетом моделирования, преследующим цель решения многомасштабных моделей и способным
использовать описание CellML, является Chaste [
FEniCS – это набор программ с открытом исходным кодом, создание которых началось в 2003 году с
целью автоматического, эффективного решения дифференциальных уравнений в частных производных [
Основа библиотеки FEniCS – это набор классов, написанных на компилируемом языке C++, которые имеют интерфейс для использования в интерпретируемом языке с неявной типизацией Python. В данном комплеске программ используется встраиваемый в синтаксис Python предметно-ориентированный язык (Embedded DSL). Например, для того чтобы решить уравнение Пуассона u = f в области c граничными условиями типа Дирихле u = u0, необходимо передать в функцию решения вариационную запись этого уравнения: WeakForm = i n n e r ( grad ( u ) , grad ( v ) ) dx
f v dx В вышеприведенном листинге переменной WeakForm присваивается следующее выражение: Z ru rvdx
Z f vdx; (3) где rg обозначает градиент функции g, а U V – скалярное произведение векторных фукнкции U и V . Равенство нулю выражения (3) для всех гладких функции с ограниченным носителем v составляет вариационную запись уравнения Пуассона. Как видно из примера, приведенный выше код во многом похож на запись формулы (3).
В библиотеке FEniCS эффективность кода на языке С++ при написании на Python удается сохранить благодаря модулю Instant. Этот модуль представляет собой динамический компилятор (just-in-time compilator), преобразующий код Python в модуль С++. Для уменьшения издержек перекомпиляции, уже готовый модуль получает свой индивидуальный номер (хеш, рассчитываемый по алгоритму SHA1). Для того чтобы иметь возможность использовать OpenMP и другие внешние библиотеки, Instant позволяет указывать заголовочные файлы, скомпилированные библиотеки, пути поиска библиотек и другие флаги компиляции.
Функции библиотеки FEniCS поддерживают параллельное выполнение как на нескольких потоках на одном вычислительном узле, так и на нескольких узлах с помощью технологии MPI. Однако, в обоих случаях требуется предварительная подготовка используемых конечноэлементных сеток.
Использование метода конечных элементов, возможность параллельного исполнения кода, хорошая документированность модулей, большое количество примеров, растущее количество библиотек, активно использующих возможности FEniCS, например, для решений задач механики сплошной среды, проблем формирования и распространения океанический волн, а также близость итогового кода к математическим обозначениям стали причиной нашего выбора этой библиотеки. 4
Вычислительные эксперименты Тестирование производилось на кластере ИММ УрО РАН "Уран". Он состоит из 178 вычислительных узлов, установленных в модулях с высокой плотностью упаковки. Вычислительные узлы оснащены процессорами Intel Xeon, работающими на частотах 2.2-3 ГГц, 16-200 ГБ оперативной памяти и графическими ускорителями NVIDIA Tesla. В общей сложности пользователям доступно 1768 вычислительных ядер CPU, 362 платы GPU и 4 ТБ оперативной памяти. Для передачи данных между вычислительными узлами используется высокоскоростная сеть Infiniband с пропускной способностью 20 Гб/сек. Для параллельного запуска на кластере установлены различные версии MPI. Дополнительно на кластер была установлена библиотека FEniCS 1.6.0. В проведенных экспериментах мы использовали MPI MVAPICH2 Intel 14.0.
Ввиду того что FEniCS поддерживает запуск в параллельном режиме, дополнительных модернизаций кода не требовалось. Для запуска вычислений с MPI использовалась команда mpiexec -f [hosts], где hosts файл с указанием вычислительных узлов и ядер, которые будут задействованы в ходе вычисления задачи.
В качестве начальных данных используется сетка, полученная путем разбиения симметричной модели левого желудочка на тетраэдры. На рис. 1 представлена схема генерации сетки. Вначале генерируются файлы полигонов в формате PLY согласно формулам, которые определяют поверхности эндо- и эпикарда. Также файлы PLY могут быть сгенерированы с помощью систем компьютерной алгебры, например, Maple. Функции из Enchanced.py удаляют повторяющиеся точки из файла PLY. Это необходимо для дальнейшей корректной работы Gmsh, который является генератором разделения расчетной области на подобластитетраэдры. Для того чтобы использовать сетку, полученную от Gmsh, в FEniCS, необходимо преобразовать ее в подходящий формат. Для этих целей предназначена утилита dolfin-convert (составная часть пакета FEniCS). Формат HDF5 можно читать в параллельном режиме, поэтому используется дополнительное преобразование из xml формата.
Для проведения тестов мы использовали сетку, имеющую 3771 вершину и 14544 тетраэдра. Программа вычисляет распространение электрической волны по левому желудочку в течение 100 безразмерных единиц времени с шагом 0.125. В начальный момент времени с помощью параметров задается зона активации.
Для оценки производительности нами была проведена серия экспериментов, в которых запуски производились на 1, 2, 4, 8, 16 и 32 ядрах. Конфигурация узлов представлена в табл. 1.
Целью эксперимента было оценить масштабируемость реализации моделирования сердца в зависимости от количества используемых ядер. Результаты представлены на рис. 2. Все тесты были запущены на одном узле, кроме теста с 32 ядрами.
В итоге, эксперимент показал хорошую масштабируемость реализации на FEniCS. Однако, ускорение замедляется с увеличением ядер из-за дополнительных затрат на передачу данных. 5 Заключение
В результате данной работы из ряда существующих программных средств моделирования свойств тел сложной формы была выбрана библиотека FEniCS. Благодаря этому наш код близок к математической Рис. 1: Схема генерации сетки Таблица 1: Конфигурация узлов кластера Параметр Значение
CPU 2 х Intel(R) Xeon(R) CPU E5-2660 0 @ 2.20GHz
RAM 94 GB
Interconnect Infiniband DDR (20 Gbit)
Operating System Linux RedHat 6.5 Рис. 2: Время моделирования в зависимости от числа использованных ядер нотации. Он транслируется в высокопроизводительный код на С/С++ с поддержкой параллельного исполнения. Помимо этого, FEniCS предоставляет широкий набор готовых инструментов для решения систем дифференциальных уравнений.
На данный момент для изучения левого желудочка сердца мы используем модель Алиева–Панфилова, состоящую из двух нелинейных уравнений в частных производных.
Реализация модели Алиева–Панфилова с помощью библиотеки FEniCS продемонстрировала масштабируемость на кластере, близкую к линейной. Тестирование производилось на кластере ¾УРАН¿ ИММ УрО РАН.
В дальнейшем мы планируем добавить поддержку таких популярных ускорителей вычислений, как GPGPU и Intel Xeon Phi, использовать более сложные и приближенные к реальности модели, которые разрабатываются в ИИФ УрО РАН, и добавить к модели механический отклик на изменение потенциала действия [12]. Благодарности
Работа поддержана программой РАН I.33П ¾Фундаментальные проблемы математического моделирования. Фундаментальные проблемы факторизационных методов в различных областях. Алгоритмы и математическое обеспечение для вычислительных систем сверхвысокой производительности¿ (проект ¾Разработка вычислительных методов и программных средств для трехмерного моделирования сердечно-сосудистой системы¿, № 0401-2015-0025). Список литературы [12] O. Solovyova, N. Vikulova, L. B. Katsnelson, V. S. Markhasin, P. Noble, A. Garny, P. Kohl, D. Noble: Mechanical interaction of heterogeneous cardiac muscle segments in silico: effects on Ca2+ handling and action potential. International Journal of Bifurcation and Chaos, 13(12):3757–3782, 2003.
Application of FEniCS framework for simulation of the cardiac left ventricle
Timofei I. Epanchintsev1;2;3, Vladimir S. Zverev1;2, Andrey V. Sozykin1;2;3 1 – Ural Federal University (Yekaterinburg, Russia) 2 – Krasovskii Institute of Mathematics and Mechanics (Yekaterinburg, Russia) 3 – Institute of Immunology and Physiology (Yekaterinburg, Russia)
Investigation of the functional characteristics of myocardium is a difficult task and requires mathematical models that describe processes on different levels. Toolkits for automated scientific computing frameworks have become popular. The aim of the paper is determine possibilities of application of such software for heart simulation on parallel computing systems. A small overview of existing frameworks (OpenFOAM, OpenCMISS, FEniCS, and Chaste) is presented. For our purposes FEniCS frameworks was chosen. We use the model of Aliev–Panfilov in order to investigate FEniCS facilities. Performance and scalability of FEniCS-based solution were studied.