научная статья по теме ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ДИНАМИЧЕСКИХ ПРОЦЕССОВ В БИОМЕХАНИКЕ СЕТОЧНО-ХАРАКТЕРИСТИЧЕСКИМ МЕТОДОМ Математика

Текст научной статьи на тему «ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ДИНАМИЧЕСКИХ ПРОЦЕССОВ В БИОМЕХАНИКЕ СЕТОЧНО-ХАРАКТЕРИСТИЧЕСКИМ МЕТОДОМ»

ЖУРНАЛ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИКИ И МАТЕМАТИЧЕСКОЙ ФИЗИКИ, 2015, том 55, № 8, с. 1380-1390

УДК 519.634

Посвящается светлой памяти А.П. Фаворского

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ДИНАМИЧЕСКИХ ПРОЦЕССОВ В БИОМЕХАНИКЕ СЕТОЧНО-ХАРАКТЕРИСТИЧЕСКИМ МЕТОДОМ1

© 2015 г. К. А. Беклемышева*, А. В. Васюков*, И. Б. Петров*, **

(*141700Долгопрудный, М.о., Институтский пер., 9., МФТИ;

**119333 Москва, ул. Губкина, 8, ИВМРАН) e-mail: amisto@yandex.ru, vasyukov@gmail.com, petrov@mipt.ru Поступила в редакцию 26.02.2015 г.

В работе представлены результаты численного моделирования механических процессов, возникающих в биологических тканях при динамических воздействиях. Для решения системы уравнений механики деформируемого твердого тела используется сеточно-характеристиче-ский численный метод на нерегулярных расчетных сетках, учитывающий характеристические свойства определяющей системы уравнений в частных производных, обеспечивающий построение более корректных вычислительных алгоритмов на поверхности раздела сред и границах областей интегрирования. Библ. 41. Фиг. 10. Табл. 2.

Ключевые слова: математическое моделирование, параллельные алгоритмы, сеточно-харак-теристический метод, биомеханика, неструктурированные сетки.

DOI: 10.7868/S0044466915080049

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

Математическое моделирование как метод исследования в данном случае обладает рядом очевидных достоинств: небольшая стоимость численного эксперимента, широта доступного диапазона изменения основных параметров, полнота получаемой в результате картины динамики скоростей, напряжений и деформаций во всем объеме рассматриваемой системы, а также области повреждения тканей (см. [1]). Эффективность этого подхода подтверждается результатами его применения в различных областях науки и техники: авиации, машиностроении, космической промышленности, сейсморазведке (см., например, [2]—[5]). Применение математического моделирования для изучения физиологических и патологических процессов, происходящих в организме человека, позволяет получить новые качественные и количественные характеристики функционирования органов в различных условиях и протекающих в них нормальных и патологических процессов. Это необходимо для прогнозирования их развития, предсказания последствий патологий, выдачи медицинских рекомендаций и разработки новых принципов диагностики на ранних стадиях различных заболеваний. Можно привести некоторые примеры успешного использования математического и численного моделирования в медицине.

Так, одними из пионерских работ по моделированию функционирования миокарда были статьи [6], [7]. Моделирование кровотока рассматривалось в работах [8]—[11], роста опухоли — в [12], процессов дыхания — в [13], [14]. Результаты расчета динамических в сосудах процессов, происходящих в глазу человека при офтальмологических операциях по удалению катаракты (фа-коэмульсификация) приведены в [15], при операциях по удалению камней в почках (литотрип-сия) — в [16]. Разработке математических моделей в иммунологии посвящена монография [17]. Исследованию распространения импульсов в волокнах Пуркинье (структуры кабельного типа,

1) Работа выполнена при финансовой поддержке гранта РНФ (14-31-00024).

проводящие импульсы в сердце) посвящены статьи [18]—[20], численному моделированию процесса залечивания ран посвящена работа [21]. Разработке реологических механико-математических моделей, описывающих поведение биологических тканей (биокомпозитов) посвящены статьи [22], [23], поведению нелинейных биологических сред — работа [24]. Применение математических методов в диагностике и анализе медицинской информации рассматривалось в работе [25], в травматологии (расчет последствий черепно-мозговых травм) — в [26]. Идеи подходов к математическому описанию психодиагностики изложены в [27], [28]. При разработке надежных математических моделей функционирования различных органов возникает ряд проблем. Наиболее сложной из них является экспериментальная верификация расчетных данных, так как соответствующие эксперименты практически отсутствуют. Более того, некоторые процессы в данной области изучать экспериментально крайне затруднительно, если вообще возможно. Нередко остаются открытыми вопросы о нахождении достоверных механических характеристик биологических тканей и о степени влияния индивидуальных различий на достоверность результатов моделирования. Решение данных проблем представляет собой актуальную задачу. Данная работа посвящена численному моделированию динамического воздействия на череп, грудную клетку и коленный сустав человека. При этом происходит распространение упругих волн в теле, и, в силу сложного внутреннего строения биологических тканей, итоговые области максимальных нагрузок формируются в результате интерференции прямых волн, а также отраженных и преломленных на внешних границах и поверхностях раздела сред.

Для математического моделирования волновых процессов в деформируемом твердом теле в приближении упругих малых деформаций и линейной упругости используется система динамических уравнений (уравнения движения и реологические соотношения) [29], [30]:

Р V; = V/Ъу + /,

1 (1)

= 1 (^Мы + ^(8**8/7 + 8п8/к))(Ук^ + Vlvk).

Здесь р — плотность среды, V — компоненты скорости смещения, а у, — компоненты тензоров напряжений и деформаций, V, — ковариантная производная по у-й координате, / — массовые силы, действующие на единицу объема, Ъу — символ Кронекера, X и ц — параметры Ламе.

Систему (1) можно переписать в матричной форме:

ди = а дц + а ди

д? х дх у ду

^ = Ах ^ + Ау д^ = / (2)

Здесь и = {vx,vy, ахх,ауу, аху}т — вектор искомых функций, х, у — независимые пространственные переменные, I — время, / — вектор правых частей, размерность которого равна размерности исходной системы, а выражения для компонент матриц Ах и Ау зависят от выбора реологических свойств исследуемой среды.

Для решения задачи (2) используется метод расщепления по направлениям. Спектральное исследование матриц Ах и Ау проведено в [31], [35], где показано, что для них существует полный набор собственных значений и собственных векторов, и приведены аналитические выражения. При переходе к инвариантам Римана система распадается на шесть одномерных уравнений. Значения инвариантов Римана переносятся с временного слоя п на временной слой п + 1 вдоль характеристических кривых (см. [2]). В данной работе используются неструктурированные треугольные сетки (см. [2], [3], [26]). В ходе расчета для каждой характеристики, выпущенной из точки на новом временном слое, по углу наклона X и шагу по времени т определяется ячейка сетки на старом временном слое, в который попала характеристика. В точке пересечения характеристики со старым временнь1 м слоем интерполируется значение соответствующего инварианта [32]. Значение инварианта переносится в точку на новом временном слое. После того, как описанным образом найдены все шесть инвариантов Римана, можно найти в ней исходные переменные и. Гибридная схема использует две опорные схемы с интерполяцией первого и второго порядка, переключаясь в зависимости от локальной гладкости численного решения (см. [33], [34]). В случае когда узел находится на границе расчетной области, рассматриваемая система уравнений имеет ровно две (см. [31], [35]) выводящие из области интегрирования характеристики, и для корректной постановки задачи требуется задание двух граничных условий для каждого граничного узла сетки: свободная граница (ат = ап = 0), заданная внешняя сила (ат = ат0, ап = ап0), заданная скорость границы (V = V,). Здесь ап и ат — нормальное и тангенциальное напряжения в

-

I

I

щг

Фиг. 1. Постановка задачи. Слева — трансверсальное сечение головы человека, справа — вид расчетной области:

1 — мозговое вещество, 2 — ликвор, 3 — внутренний слой компактной костной ткани, 4 — слой губчатой костной

ткани, 5 — внешний слой компактной костной ткани.

граничной точке, сти0 и стт0, и0 — заданные извне нормальное и тангенциальное напряжения и вектор скорости в точке границы. Расчет контактной границы между двумя телами в целом аналогичен расчету границы тела. Однако в точке контакта присутствуют два узла — по одному из каждого контактирующего тела, и уравнения связей должны задавать 4 условия. В данной работе используется ряд контактных условий: скольжение тел друг относительно друга (Vп = Vп,стп = стп,стт = стт = 0) и слипание тел (^п = = V). Здесь а„ и стт — нормальное и тангенциальное напряжения в граничной точке. Символы с волной относятся к первому телу, без волны — ко второму.

Из нейрохирургической практики известно, что области поражения мозга при черепно-мозговой травме не всегда совпадают с областями, прилежащими к месту удара (см. [26], [36]). По результатам анализа ряда работ можно сделать вывод о том, что большинство исследовательских работ в области моделирования ЧМТ было проведено за рубежом с применением различных вариаций методов конечных элементов (МКЭ); обзор работ, посвященных численному исследованию ЧМТ, выполнен в [1], [26]. Несмотря на то что моделирование головного отдела человека с помощью МКЭ в последние десять лет интенсивно развивается, оно еще далеко от того, чтобы иметь возможность объяснить механизмы повреждения мозга и предсказать последствия ЧМТ. Было проведено некоторое количество сравнительных исследований с использованием различных контактных условий (см. [37]). Все они приводят к заключению, что эффект ударного воздействия на голову человека чувствителен к способу моделирования контактного взаимодействия череп-мозг. Сеточно-характеристические численные методы, ориент

Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.

Показать целиком