Автоматизированный подход к параметрической генерации модели и МКЭ-оценке процесса ее сжатия: пример каркаса транскатетерного протеза клапана сердца

С распространением приобретенных клапанных     высокотехнологичных

малоинвазивных

заболеваний развиваются методы их лечения: от транскатетерных имплантаций с 2007 года [2]. К открытых операций с 1952 года [1] до     преимуществам последних относят

значительно

0000-0003-2404-2873

0000-0003-3211-1250

0000-0001-6099-0315

0000-0001-7477-3979

Эта статья доступна в соответствии с условиями лицензии Creative Commons Attribution-NonCommercial 4.0 International

License (CC BY-NC 4.0)

This work is licensed under a Creative Commons Attribution-NonCommercial 4.0 International License (CC BY-NC 4.0)

меньшие риски кровотечений, инсультов, более быстрое восстановление, а также выгода с точки зрения экономических затрат [3,4]. Именно поэтому в России и в мире транскатетерное протезирование переживает расцвет. На территории РФ с 2018 по 2023 год было зафиксировано увеличение в 1,91 раза количества таких вмешательств: с 920 из 10673 [5] до 1756 из 11404 [6] случаев в год. Растущая потребность в заместителях такого типа и развитие численных подходов [7–9] требует от инженерных и научных коллективов разработки собственных подобных конструкций, в т.ч. для снижения зависимости от импорта.

Конструктивно, транскатетерные протезы клапанов сердца содержат три основных компонента:

• 1.    Ячеистый каркас из сплавов никелида титана [10], кобальт-хрома или нержавеющей стали;

• 2.    Створчатый аппарат, выполненный из перикардиальной ткани животных или полимерных соединений [11];

• 3.    Уплотняющая юбка-манжета, снижающая утечку крови между протезом и окружающими тканями.

Помимо самого транскатетерного протеза для имплантации требуются устройство сжатия – кримпер – и длинная (не менее 1 метра) гибкая система доставки. Процесс имплантации в аортальную позицию (transcatheter aortic valve replacement – TAVR) происходит с доступом через крупные кровеносные сосуды – бедренную, реже подключичную или сонную артерии [12; 13], без открытого доступа к сердцу. Контролируя процесс доставки с использованием методов визуализации в реальном времени (флюороскопия), протез позиционируют в пораженном клапане сердца и расширяют изнутри баллоном до необходимого диаметра или снимают ограничивающий кожух для самораскрывающихся протезов. Перед проведением такой процедуры протез необходимо разместить на системе доставки так, чтобы он обладал достаточно малым диаметром 5-7 мм для прохождения по сосудам пациента, причем недостаточное сжатие приведет к сопутствующим нежелательным повреждениям окружающих тканей во время доставки к сердцу.

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

• 1.    Сжимаемость опорного каркаса без разрушения, значительных искажений геометрии.

• 2.    Возможность обратно восстановить свою форму в области имплантации – пораженном клапане сердца, без разрушения.

• 3.    Обеспечение надежности фиксации в зоне имплантации – каркас должен иметь достаточную радиальную жесткость, чтобы избежать деформации окружающими тканями и при этом не сместиться из просвета клапана сердца, т.к. процедура не предполагает фиксации протеза швами [14; 15];

• 4.    Обеспечение долговечности каркаса после цепочки «сжатие – расширение баллоном» с сохранением геометрии, предотвращая появление зон концентрации напряжений, являющихся потенциальными очагами разрушения конструкции в условиях долгосрочной циклической нагрузки. При этом только за 5 лет работы в организме протез подвергается 200 млн циклов нагрузки [16].

Ввиду исключительной важности и конструктивной сложности опорного каркаса в составе транскатетерного протеза, научные коллективы ведут активный поиск методов проектирования его «оптимальной» геометрии, в том числе с использованием аддитивных технологий [17]. Из всех них можно выделить последовательность шагов «генерация дизайна» – «численное моделирование» методом конечных элементов (МКЭ) [18–21]. Она может быть дополнительно модифицирована алгоритмами оптимизации (например, NSGA-II [18]), или использовать экспертные оценки [22] для выбора наилучшей геометрии опорного каркаса. Основные работы в рамках данной задачи улучшения дизайна исследуемого объекта подходят со стороны варьирования нескольких параметров в рамках единого дизайна [18; 20; 22–25].

Ранняя из указанных выше работ в данной области датирована 2013 годом, когда Li H. [22] с соавторами произвели оптимизацию баллонорасширяемого коронарного стента. Оптимизационная процедура была основана на варьировании 5 геометрических параметров, описывающих толщины перемычек при фиксированной длине и диаметре, а также расстоянии между прямоугольными ячейками с целью оптимизации дизайна стента для его безопасного раскрытия в зоне атеросклеротического поражения сосуда. В 2021 году за авторством Russ J.B . [23] было опубликовано исследование, нацеленное на оптимизацию баллонорасширяемых стентов для протезосодержащих кондуитов, применяемых в детской кардиохирургии. Она нацелена на получение дизайна путем варьирования трех параметров при фиксированной длине и диаметре: толщина и ширина струта и количество ячеек. Форму ячейки можно описать как область между двумя синусоидальными линиями, отличные знаком, однако дизайн каркаса имеет всего один ряд (при длине 20 мм). Здесь авторы с применением NSGA -II улучшали качество расширения баллоном стента, оценивая диаметры сосуда в области имплантации, пиковые радиальные силы, метрику Гудмана и гидростатический критерий для прогнозирования разрушения. В работе S. Barati (2022 год) [18] описывают переложение таких методик на сложные модели геометрий каркасов TAVR протезов. Авторы используют 6 параметров дизайна, полностью состоящего из ячеек равного размера: их количество по высоте и диаметру (форма аналогична Russ J.B. [23]), толщину и ширину струта, проксимальный диаметр и число, характеризующее искривление формы. Описанный в этой работе подход ориентирован на оптимизацию самораскрывающегося

Рис. 1. Блок-схема работы алгоритма

каркаса, выполненного из никелида титана с применением NSGA -II для четырех целевых функций качества функционирования: максимальная сила сжатия, радиальная сила, площадь контакта с окружающими тканями (анкировка) и диаметральные искажения после установки в область имплантации. Схожая работа была проведена Kovarovic B. (2022) [20] в целом описывающая опыт оптимизации всех компонент полимерного протеза TAVR . В отличие от исследования S. Barati , каркас содержит только 3 ряда базовых ячеек ромбовидной формы с вводной стороны, занимающие примерно 30% высоты всего дизайна. К сожалению, авторы не раскрывают подробностей о варьируемых параметрах и оптимизационном процессе, но говорят о том, что уменьшение толщины и ширины струтов позволило значительно улучшить фиксирующую способность всего каркаса без потерь в возникающих эффектах напряженно-деформированного состояния.

С нашей точки зрения несмотря на то, что вышеописанные работы показывают актуальность данной тематики, основным ограничением является малое количество геометрических параметров построения дизайна. Для преодоления этого в данной работе будет представлен подход с расширенным набором параметров для генерации параметризованной модели каркаса транскатетерного протеза клапана сердца с использованием пакета CadQuery [26], обеспечивающий генерацию геометрии с последующей интеграцией в среду МКЭ-анализа Abaqus/CAE 2021 (Dassault Systemes, Франция) для проведения симуляции процесса его сжатия.

Материалы и методы

Учитывая сложность предъявляемых требований к функционированию опорного каркаса транскатетерного протеза в настоящей версии алгоритма, мы ограничили его функционал первой и самой критической нагрузкой - имитацией его сжатия до диаметра 6 мм с оценкой напряженно-деформированного состояния.

Разработанный подход реализован на языке программирования Python 3.11 и состоит из трех ключевых блоков: «Генератор», «Численный анализ», «Оценка» (рис. 1).

Все параметры для запуска считывали из конфигурационного файла формата yaml и получали словари, содержащие информацию о:

• •    именах варьируемых параметров геометрии;

• •    значениях и диапазонах геометрических параметров геометрии;

• •    модели материала, применяемой в блоке «Численный анализ»;

• •    настройке решателя (например, количество ядер для запуска расчета, пути расположений, список полей для экспорта из результатов оценки напряженно-деформированного состояния и т.д.).

Всего для описания одного дизайна применяли 11

Рис. 2. Принцип генерации геометрии каркаса: а – элементарная ячейка; б – размещение ячеек разных уровней по высоте

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

• 1.    Определить перечень ограничений параметров, например, по диаметру каркаса.

• 2.    Определить диапазоны варьирования оставшихся параметров.

• 3.    Описать модель материала, из которой будет выполнен предполагаемый каркас.

• 4.    Указать настройки решателя.

• 5.    Запустить.

• 6.    Оценить полученные результаты.

Далее опишем основные идеи и их реализацию, заложенные в каждый из трех блоков последовательности «Генератор» – «Численный анализ» – «Оценка». Далее будем применять следующую нумерацию параметров: число в верхнем индексе геометрического параметра означает номер уровня базовой ячейки, иначе он относится ко всему дизайну.

Данный блок отвечает за создание дизайна оцениваемого каркаса и принимает на вход словарь геометрических параметров и имя итогового файла дизайна. Функция создания использует библиотеку CadQuery версии 2.4.0. В ее основе лежит разработанный универсальный для данной задачи базовый генератор (рис. 2, а ). С его использованием производили создание ячеек трех уровней, отличающихся пропорциями, для создания сетчатого «рисунка» каркаса (рис. 2, б ) и распределяли по окружности, учитывая количество их повторов N Repeat . Затем, из цилиндра определенных диаметра (реализация позволяет как использовать его в качестве параметра построения, так и оставаться постоянным значением) и толщины стенки – будущего струта – вычитали

Таблица 1

Диапазоны варьиования и фиксированные значения геометрические параметров построения дизайнов

Название Описание Диапазон / значение D Диаметр каркаса 26 мм h1 Длина нижнего прямого участка [0,1; 1] мм h2 Длина среднего прямого участка [0,1; 1] мм h3 Длина верхнего прямого участка [0,1; 1] мм h3 Длина среднего прямого участка ячейки 3 уровня [1; 6] мм width Ширина нижнего участка базовой элементарной ячейки 0,5 мм hcell Общая высота ячейки 1 уровня [1; 10] мм NRepeat Количество повторов шаблона распределения ячеек по окружности 12 А Коэффициент ассиметрии ячейки [0; 1] pad Ширина струта (расстояние между ячейками) 0,5 мм THK Толщина стенки каркаса 0,5 мм полученные ячейки. В рамках данной работы было принято, что будут рассмотрены только дизайны каркасов с NRepeat = 12 повторов, однако, реализация подхода позволяет использовать любые целые значения больше 3.

Базовый генератор ячейки принимает набор геометрических параметров (рис. 2А) и, затем, начинает построение геометрии обходом по часовой стрелке с

Рис. 3. Схема блока «Численный анализ» очередности действий создания численной сборки сжатия каркаса

левой нижней точки. Значения длин l1 и l2 являются вычисляемыми по следующим формулам, полученным эмпирически, исходя из логики, что если коэффициент асимметрии A = 1, то 1 1 = 12:

1 1 =

hi ell - ( ^h 1 + ^h 2 + h 3 )

A ,

• ¹ 2 = ^h Cell - ( h l + ^{h 2} + ^h 3 ) - ^l 1 ,

где h ell - соответствующее значение высоты ячейки 1 уровня (рис. 2А); h 1 - h ₃ - высоты прямых участков (рис. 2 А); A - коэффициент отличия длин 1 1 и 1 2. В рамках одного дизайна ячейки 2 и 3 уровня отличались длиной участка h 2 с сохранением всех остальных параметров неизменными.

Полученные на данном этапе ячейки 1–3 уровней распределяли по вертикали для образования «шахматного» порядка (рис. 2Б), а затем вычитали из оболочки цилиндра определенного диаметра с толщиной стенки, имитируя технологический процесс производства - лазерную резку из трубки [27,28]. В результате работы данной функции получали трехмерный дизайн каркаса с последующим сохранением в файлы геометрий формата STEP.

В рамках данной работы для качественной оценки состоятельности подхода использовали следующие диапазоны и значения параметров (табл. 1):

На данном этапе производили создание сборки численного анализа.

Вызванный скрипт в среде Abaqus/CAE 2021 выполняет протокол создания численной сборки, схожий с порядком действий оператора программы в интерактивной оболочке Abaqus/CAE (рис. 3). На данном этапе сначала импортировали модель каркаса, создавали участок цилиндра, имитирующего сжимающее устройство, который отстоит на 0,2 мм от исследуемой модели. В качестве базового элемента сетки были использованы шестигранники C3D8. Сетка цилиндра состояла из поверхностных прямоугольных элементов типа SFM3D4. Поскольку в рамках данной работы было принято NRepeat =  12, то для ускорения процесса вычисления использовали симметричную постановку, включающую 1/NRepeat = 1/12 от всего каркаса, причем медианные значения получаемых расчетных сеток составили 10130 (квартили 25%: 9041; 75%: 11345) элементов для всей постановки. В связи с этим необходимо было зафиксировать грани «рассечения» (u2=0 в цилиндрической системе координат). Расчет проводили на промежутке 0-1 секунда, где сжатие происходило на отрезке 0-0,75 секунд, а остаток времени – оценка «рекойла» каркаса.

В качестве материала каркаса задавали нержавеющую сталь 316 с модулем Юнга 233000 МПа, коэффициентом Пуассона 0,3 и плотностью 7,98∙103 кг/м3. Пластическую компоненту модели задавали с использованием данных от производителя: 414 МПа для наступления фазы пластичности и предел прочности 933 МПа при 45% удлинения.

В данной постановке были использованы два вида взаимодействия компонентов сборки: (а) сектора цилиндра и внешняя поверхность каркаса; (б) граней каркаса между собой (« self-contact »). В обоих случаях использовали контакт твердых тел (« Hard contact ») с коэффициентом трения 0,2.

Результатом на данном этапе становились входные inp-файл сборки, содержащие все необходимые данные для запуска численного расчета. Затем, их передавали Abaqus для численного анализа солвером Standard . После этого, полученные файлы результатов симуляции сжатия в формате odb передавали в последующий блок «Оценка».

Данный этап автоматически выполняли после завершения расчета процесса сжатия. С использованием скрипта на языке Python 2.7 и пакетов API Abaqus/CAE производили расшифровку и сохранение данных расчета в формате csv в моменты времени 0,75, 1 секунду или в последний сохраненный момент (в случае, если возникли ошибки при численном анализе). В данном блоке оценивали три параметра: максимальное напряжение по Мизесу ( S , МПа), сумму сил реакции опоры ( RF , Н), которую брали как первую компоненту тензора сил в цилиндрических координатах (ось Z соосна каркасу) в моменты времени T = 0,75 сек, 1 сек и последний стабильный (для случаев «развала» решения или превышения временного лимита 1 час на дизайн). Рекойл каркаса после снятия внешнего воздействия ( Rc , мм) вычисляли как разность диаметров. Результаты записывали в таблицу, хранящая в себе данные о геометрических параметрах построения, результатах симуляции сжатия, а также затраченное время на один дизайн.

С точки зрения значений этих параметров их можно соотнести с прикладной интерпретацией для предполагаемых пользователей:

• 1.    Максимальное напряжение по Мизесу с запасом прочности материала вследствие возникающих динамических эффектов.

• 2.    Сумма сил реакции опоры с тем насколько эффективно имплантированный клапан будет сопротивляться сдавливанию со стороны окружающих тканей сердца и сохранять свою форму после того, как баллон сдуют. Так, высокие значения требуют приложения значительных усилий для сжатия и, как следствие, ведут к невозможности монтажа в зоне имплантации баллонами, обеспечивающими давление 5-7 атмосфер [28], а низкие – могут не обеспечить достаточную силу фиксации.

• 3.    Рекойл каркаса показывает эластическую отдачу материала, который может затруднить как размещение протеза на системе доставки ввиду недостаточной фиксации, что увеличивает риски травмирования сосудов во время доставки.

Дополнительно вычисляли геометрические параметры, такие как длина получаемого каркаса ( L , мм), поскольку в явном виде мы не используем этот параметр и его можно применить как один из критериев ограничения – constrain – оптимизационного процесса в будущем.

Весь алгоритм был запущен на ноутбуке Asus UX3402V ( Intel Core i5-13500H , 16 Гб ОЗУ) под управлением операционной системы Ubuntu 24.04.03 с режимом энергопотребления «Производительный».

Оценка результата работы алгоритма

В рамках данной работы был использован метод оптимизации Random Search , заключающийся в случайном поиске наилучшего решения. Такой подход не является оптимальным для решения практических задач и нацелен на оценку возможностей описываемого подхода в целом ввиду покрытия большого пространства возможных комбинаций значений параметров построения. Всего в его рамках было получено 1225 уникальных дизайнов каркасов транскатетерных протезов клапанов сердца. Валидный набор полей « S , RF » присутствовал у 18/450 (4,7%) на момент T =0,75 c и у 36/1225 (2,9%) на момент T =1,0 c, медианное (квартиль1; квартиль3) затраченное время на генерацию и анализ одного дизайна составило 3 мин 40 сек (2 мин 40 сек; 5 мин 11 сек). Такая низкая доля валидных решений указывает на ограниченную численную устойчивость текущей постановки, являющийся следствием автоматизированной генерации расчетной сетки. Показанные далее примеры «удачных» дизайнов не следует трактовать как исчерпывающие данные, а как скрининговые.

Приведем некоторые результаты численной оценки каркасов, способных быть сжатыми до диаметра 6 мм (рис. 4). На изображении приведены данные пяти дизайнов, которые по результатам численного анализа способны быть сжатыми до конечного диаметра.

Рис. 4. Результаты и внешний вид некоторых дизайнов, способных быть сжатыми до диаметра 6 мм

Рис. 5. Результаты анализа дизайнов 229 и 750. А: Визуализация результатов численного анализа и рекойл, после снятия внешнего воздействия; Б: график зависимости необходимой для сжатия силы RF от диаметра каркаса

Рис. 6. Демонстрация возможности разработанного подхода создавать дизайны каркасов различной формы

Важно отметить, что при схожих значениях напряжения по Мизесу (933,00 и 879,20 МПа) и длины каркаса (25,40 и 24,08 мм) для дизайнов 229 и 750 (рис. 5А) наблюдали различие в RF – пиковые значения 1374,20 Н и 284,24 Н (рис. 5Б), обусловленное контактом стенок ячейки для первого варианта (рис. 5А). Как известно из литературных данных [29], сила, необходимая для сжатия каркаса, варьируется в диапазоне от 50 до 350 Н. Полученные экспериментальные данные свидетельствуют о том, что вариант геометрии каркаса ID 229 будет сложно подвергнуть сжатию до 6 мм и возникающие напряжения равны пределу прочности материала, что делает его несостоятельным и при отборе/оптимизации должен быть исключен из рассмотрения.

Далее продемонстрируем технические способности разработанного подхода генерации различных дизайнов каркасов транскатетерных протезов клапанов сердца. Так был создан набор случайных геометрий (рис. 6), диапазоны варьирования параметро которых были выбраны следующим образом (табл. 2).

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

Ограничения

Несмотря на полученные результаты, данное исследование имеет ряд ограничений, которые необходимо учитывать при интерпретации данных. Прежде всего, использование симметричной постановки задачи с моделированием 1/12 части конструкции вносит определенную погрешность в сравнении с расчетом полномасштабной модели каркаса. Кроме того, предложенный алгоритм был ограничен использованием одного базового типа элементарной ячейки. Хотя варьирование геометрических параметров позволило получить широкое разнообразие конфигураций, архитектурные решения, основанные на принципиально иных топологиях или ячейках сложной формы, остались за рамками текущего анализа.

Механическое поведение конструкций оценивалось в условиях симметричного сжатия сектором цилиндра, что не в полной мере отражает специфику работы реальных систем доставки. Устройства радиального сжатия, применяемые в инженерной и клинической практике, воздействуют на каркас плоскопараллельными площадками, формирующими многогранник, из-за чего распределение локальных нагрузок может отличаться от модели с цилиндрическим воздействием. Наконец, работа была сфокусирована исключительно на этапе сжатия каркаса, в то время как вопросы долгосрочной усталостной прочности и взаимодействия протеза с тканями корня аорты в условиях динамического кровотока требуют дальнейшего углубленного изучения

Заключение

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

Таблица 2

Диапазоны варьиования и фиксированные значения геометрические параметров построения дизайнов для демонстрации возможностей

D, мм h1, мм h2, мм h3, мм h23, мм width, мм 26 [0,1; 2] [0,1; 2] [0,1; 2] [1; 6] 0,5 hcell , мм Nrepeat, разы А pad, мм THK, мм [1; 20] [10;20] [0; 2] 0,5 0,5 диаметра 6 мм (~18 Fr). Реализован и описан способ связи создания трехмерных моделей в CAD-среде на языке программирования Python 3.11 и его интеграция с Abaqus / CAE  через открытые API.  Важнейшим критерием новизны описанного подхода является возможность разработки новых видов опорных каркасов транскатетерных    протезов    клапанов    сердца, ограниченного    только    формой    заложенной элементарной ячейки. Разработанный алгоритм предоставляет инструмент изучения большого пространства    геометрических    параметров    с автоматической        оценкой        напряженно - деформированного состояния при моделировании сжатия. Полученные результаты работы описанного подхода показывают различие в возникающих эффектах напряженно-деформированного состояния при сжатии. Такой подход «генератор»-«численный анализ»-«оценка» является основой для построения оптимизационного    процесса,     например,     с использованием многокритериального Байесовского оптимизатора (MOBO - Multi-Objective Bayesian Optimizer) или генетического алгоритма.

Актуальная версия доступна для ознакомления по ссылке:

Приложение

Основным ограничением прямого использования API Abaqus/CAE - несовместимость версий Python (2.7 в Abaqus/CAE и Python 3.11 в данной работе). Поэтому было принято решение о создании копии конфигурационного файла yaml в формате json, который автоматически перезаписывали на каждой итерации алгоритма. Для этого были использованы стандартные пакеты os, json, tempfile. Код запуска создания сборки численного анализа реализован на языке Python2.7   с использованием пакетов   numpy, math и пакетов API abaqus. Из-за несовместимости версий, последний скрипт вызывали с использованием функции subprocess.run (cmd, shell=True), которая передает в консоль операционной системы строку команды:

cmd = 'abaqus cae noGUI="{script_path}" --"{json_path}"',