Моделирование спиральных волн в аорте Текст научной статьи по специальности «Физика»

УДК 531/534:[57+61]

МОДЕЛИРОВАНИЕ СПИРАЛЬНЫХ ВОЛН В АОРТЕ

В.А. Батищев1, Д.С. Петровская1, Ю.А. Устинов2

1 Кафедра теоретической и компьютерной гидродинамики Южного федерального университета, Россия, 344006, Ростов-на-Дону, ул. Мильчакова, 8а, e-mail: [email protected]

2 Кафедра теории упругости Южного федерального университета, Россия, 344006, Ростов-на-Дону, ул. Мильчакова, 8а, e-mail: [email protected]

Аннотация. В работе рассчитаны два типа спиральных волн, распространяющихся в аорте - короткие спиральные волны и квазистационарные моды. Исследование проведено на основе системы Навье-Стокса и динамических уравнений тонкой упругой оболочки. Аорта моделируется круговым цилиндром, ограниченным изотропной оболочкой. Показано, что короткие спиральные волны заполняют все поперечное сечение цилиндра в отличие от длинных спиральных волн, которые локализованы в пограничном слое вблизи оболочки. Упругие свойства оболочки оказывают слабое влияние на короткие волны. Часть спиральных мод локализуется в критическом слое вблизи оси цилиндра. Механизмом переноса как коротких спиральных волн, так и квазистационарных мод является стационарный поток.

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

Введение

Во второй половине прошлого века появились сообщения об обнаружении винтовых течений крови в артериях человека и животных [1, 3, 6, 7, 13, 14]. В работе [7] на основе анализа экспериментальных данных показано, что закрученный поток крови формируется в левом желудочке сердца. К механизмам его возникновения относятся «ротационная деформация миокарда левого желудочка», сокращение трабекулярных мышц желудочка, придающих «потоку изначальное вращение» [12, 13], а также «аппарат митрального клапана, папиллярные мышцы и маятникообразное (качательное) движение сердца в систолу» [5, 10]. Закрученное течение крови в сердечном желудочке порождает спиральное течение в восходящей части аорты. В [7] отмечено, что направление вращения крови в аорте «совпадает с направлением вращения крови в левом желудочке». Однако в некоторых случаях возникает обратное течение, что вызвано «индивидуальными морфофункциональными особенностями выходного отдела левого желудочка». Итак, в этих работах на основе анализа экспериментальных данных показано, что одной из причин возникновения спиральных течений крови может быть закрученная структура стенок левого желудочка сердца, что порождает вихревые течения крови на входе в аорту и перенос вихрей потоком жидкости. Среди причин возникновения винтовых течений крови в аорте могут быть и механические свойства стенок кровеносных сосудов. В [2, 11] показано, что винтовая

© Батищев В.А., Петровская Д. С., Устинов Ю.А., 2013

Батищев Владимир Андреевич, д.ф.-м.н., профессор кафедры теоретической и компьютерной гидродинамики Южного федерального университета, Ростов-на-Дону

Петровская Дарья Сергеевна, магистрант кафедры теоретической и компьютерной гидродинамики Южного федерального университета, Ростов-на-Дону

Устинов Юрий Анатольевич, д.ф.-м.н., профессор кафедры теории упругости Южного федерального университета, Ростов-на-Дону

анизотропия стенок сосудов мышечного типа [9] приводит к возникновению длинных спиральных волн, причем эти волны локализуются в пограничном слое вблизи стенок кровеносных сосудов. В данной работе рассчитаны как короткие спиральные волны, так и квазистационарные моды, заполняющие все поперечное сечение цилиндрического сосуда, моделирующего аорту. Эти волны слабо зависят от упругих свойств стенок сосудов. Механизмом переноса коротких волн является стационарный поток, который наблюдается в аорте. В [8] приведены значения скорости этого потока для аорты собаки. Длинные и короткие спиральные моды изменяют направление вращения жидкости на обратное либо по сечению, либо со временем в течение сердечного цикла. В статье рассчитаны также и квазистационарные спиральные моды, в главном приближении не зависящие от времени, причем часть мод не изменяет направления вращения жидкости.

Постановка задачи

Рассматривается задача о движении спиральных волн малой амплитуды в аорте, в которой распространяется стационарный поток жидкости исчезающей вязкости. Аорта моделируется круговым цилиндром, ограниченным тонкой упругой изотропной оболочкой толщины h и радиусом срединной поверхности, равным а. Спиральные волны распространяются на фоне стационарного потока и длинных продольных пульсовых волн, которые наблюдаются в крупных кровеносных сосудах [8]. Движение жидкости описывается системой уравнений Навье-Стокса

д v

——+ (v,Vv) = -p-1Vp + иAv, divv = 0.

Здесь v - вектор скорости, v = (vr,ve, vz ) ; (r,9,z) - цилиндрические

координаты; p - плотность; p - давление; и - кинематический коэффициент вязкости. Течение жидкости считается периодическим по времени с заданным периодом T. Значения параметра T приведены в [8]. Рассматривается осесимметричная задача, т.е. поле скоростей, давление и смещение точек срединной поверхности оболочки не зависят от азимутальной координаты 9. Уравнения движения приводятся к безразмерному виду, в котором в качестве масштабов длины, времени, скорости, давления и смещений точек оболочки приняты соответственно следующие параметры:

а, —, U, pU2, U. Здесь ю = —, U - максимальное значение скорости частиц ю ю T

жидкости в поле длинных продольных волн, p - плотность жидкости. Длинные

продольные пульсовые волны, распространяющиеся в жидкости, заполняющей

цилиндр, ограниченный тонкой упругой оболочкой, изучались во многих работах [2, 4,

8, 9, 11]. В [2, 11] исследованы волны с учетом винтовой анизотропии оболочки,

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

которая описывает стационарное движение жидкости и моделируется течением

Пуазейля. В связи с этим предположим, что стационарный поток определяется полем

скоростей, у которого осевая компонента скорости изменяется квадратичным образом

по радиальной координате

R

vz = V0(r ) = Rp (1 - r 2)/R vr = v9 = 0 P = c - 4^ z^7- (1)

R

U Up

Здесь R =--------, Rp =——, где U p - характерная скорость стационарного потока,

(юа) p (roa) p F

2Q ^

U = —-. Q - среднии по времени расход жидкости в цилиндре. sv - малый параметр, па

^v =

v 1

, —т- (для аорты собаки еу « — ). Отметим, что функции (1) есть точное решение Шгоа2) 13

системы Навье-Стокса. Окружная компонента скорости удовлетворяет уравнению

д v9 д t

+ R

д v9 д v9 vrv9

v —в- + v —- +

д r д z r

= є2| Vч - Ц-1. (2)

у

r

Динамические уравнения тонкой изотропной оболочки в рамках безмоментной теории для окружной компоненты смещения с учетом гидродинамического воздействия на срединную поверхность оболочки (г = 1) приводятся к безразмерному виду

'2 Р h д2 и ( д V Ц ^

д u 2 Ро h , ч д и

- є

•л 2 k \ о / Л .

д z р а д t

(1+ Vо= є2 є2 (1 + Vо)

9 9

д и0 2а р

ц0 = —-, 8к = гоа-----.

0 д г к Чье

Здесь и0 - окружная компонента вектора смещений точек срединной поверхности оболочки; Е, р0 - модуль Юнга и плотность материала оболочки соответственно; у0 - коэффициент Пуассона, который принимается равным 0,5. Плотность р0 считается близкой к плотности жидкости. Параметр 8к мал (для аорты собаки 8к ~10-2 [8]).

Асимптотические разложения

Решение задачи представим в виде суммы двух векторов V (г, г, г1, г) = У1 + У2. Здесь г1 - медленная осевая координата, г1 = 8кг . Вектор У1 = (ц , ц0 , , р, иг, и0, иг)

описывает стационарное течение (1) и длинные волны, зависящие только от координат (г1, г, г). Вектор У2 = (, (0, ц, q, и , и0, иг ) описывает короткие волны малой

амплитуды. Первые четыре координаты векторов У1 и У2 - это компоненты векторов

скорости и давления, а последние три - компоненты смещений точек оболочки.

Задача о распространении длинных продольных волн в жидкости, заполняющей цилиндр с упругой стенкой, изучалась для различных случаев во многих работах [2, 4, 8, 9, 11]. Асимптотические разложения компонент вектора У1 представим в виде

рядов по степеням малого параметра 8к :

= Цг0 + У (Г) + 8кЦг1 + Ц0 = 8кЦ00 + ... (8к ^0). (3)

Остальные компоненты разлагаются в аналогичные ряды. Здесь предполагается, что амплитуда окружной компоненты скорости мала и имеет порядок ц0 = 0(8 к). Главные

члены асимптотических рядов для координат вектора У1 известны, а их значения для 1ББМ 1812-5123. Российский журнал биомеханики. 2013. Т. 17, № 1 (59): 55-63 57

изотропной оболочки приведены, например, в работах [4, 8]. В [2, 11] изучены длинные волны с учетом анизотропных свойств оболочки. Длинные спиральные волны, описываемые компонентой ц0 , локализуются в пограничном слое вблизи поверхности

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

Короткие спиральные волны

В отличие от длинных спиральных волн, короткие спиральные волны заполняют все поперечное сечение цилиндра и движутся со скоростью стационарного потока. Эти волны описываются компонентой ^ вектора У2. Предположим, что ^ имеет такой

же порядок малости 0(8к ), что и длинноволновая компонента ц0 . Порядки остальных

компонент вектора У2 находятся из уравнений Навье-Стокса и динамических

уравнений оболочки. Асимптотические разложения компонент вектора У2 представим

в виде

^0 = 8Л + ... ,

^ =^к2+ ..., (8к ^ 0).

и0 = 8 к и01 + ....

Аналогичные разложения записываем и для остальных компонент вектора У2.

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

приближении компонента w01 удовлетворяет на стенке условию прилипания ^01 = 0 (г = 1), а упругие свойства оболочки проявляются только в высших

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

+ к ((Г, Г„ 1) + <>)))- = ¿2 - ^0,

д 1 д 2 ^ г ) (4)

w0l = 0(г =1), w0l = 0 (г=°Х w0l ^0(2 ^да).

Здесь медленная переменная г1 входит в (4) как параметр. По времени

выполняется условие периодичности.

Решение задачи (4) построим в виде отдельных мод. Учитывая, что

коэффициенты в (4) зависят только от медленной осевой координаты г1 и не зависят от

2 , разделяем переменные w01 = Ж01 (г, 1, г1) в'ь. Для функции Ж01 сформулируем

краевую задачу, которая получается из (4) заменой производных по координате 2 на множитель ¡к.

Длинная продольная волна при малых 8У имеет характер пограничного слоя вблизи оболочки, поэтому функция 0 представима в виде 0 = ц°°(1 , г1) + 0( 51,1, г1),

где 51 = (1 - г)/ 8У . Функция ц°° описывает решение вне области пограничного слоя, а Ь2° - поправка к ц°°, устраняющая невязку в граничном условии при г = 1. Функция Ь° имеет характер пограничного слоя вблизи оболочки. Функция ц°° не зависит от

радиальнои координаты r , а ее среднее значение по времени равно нулю, поэтому функцию W01 представим в виде

W = Q(t,zi)W, Q = exP(-ikRjV0O(t,zi)dt)• (5)

Теперь, применяя метод пограничного слоя, решение строим в виде суммы

W = w + H, где w - значение функции W вне области пограничного слоя, а

H - погранслоИная поправка к w. Учитывая, что вне области пограничного слоя

функция hz0 = O, приходим к уравнению для функции w :

+ ikRV0(r) w = &V

д t

(02w 1 дw ( 2 1 ^ ^

—7 +-----------Ik2 + — I w

д r r д r I r J

Краевые условия для w используем те же, что и в (4). Далее разделим переменные w = F(r)exp(-int). Функция F(r) и комплексное волновое число к определяются из краевой задачи на собственные значения

d2F 1 dF ( 2 1 Л zr i ( , г> л 2чЛт7

+-| к2 + — IF = — I -n + kR (1 - r2) IF

dгí г dг ^ г2) 82 ^ р ') ' (6)

F (0) = 0, F (1) = 0.

Результаты расчетов

Приведем выражение для спиральных мод:

^1 = Опт(1,21)((г)ехР(((к( - п1)) + Нп т ехрО!^)). (7)

Здесь п,т = 1,2,3,.... Функция Fnm и параметр кпт определяются из краевой задачи (6). Функция Опт находится по формуле (5) , в которой к следует заменить на кп т. Первое слагаемое в (7) описывает решение вне области пограничного слоя. Функция Нп т - погранслойная поправка, локализованная вблизи оболочки. Задача (6) решалась численно методом пристрелки, причем вблизи оси цилиндра проводилось сращивание со степенной асимптотикой. Комплексное волновое число кпт = кг + iki

зависит от параметров Яр, 8Ц, п. Численные расчеты показали, что логарифмический декремент затухания к{ и волновое число кг монотонно убывают с ростом Яр. Для

т = 1, п = 5, Я = 2 и 8 =—1— приведем значение к = 2,7347 + i0,2872. Введем р " 13,116 11

обозначение Fnm = gnm + ¡/пт. Расчеты показали, что при т = 1, п >1 с ростом индекса п увеличиваются значения параметров кг и к, а функции gnm и /пт локализуются вблизи оси цилиндра, причем локализация происходит и при 8У^ 0. На рис. 1 изображена зависимость функций g51 (кривая 7) и /51 (кривая 2) от радиальной координаты г при п = 5, т = 1. Физические параметры, выбранные для расчетов, соответствуют кровеносным сосудам собаки [13].

0,4

0

0,75

г

1

Рис. 1. Спиральная мода

При ^ 0 вблизи оси цилиндра возникает вязкий критический слой для

т = 1, п > 1 с толщиной порядка О ((7) • При выходе из этого слоя спиральные моды

затухают по экспоненциальному закону. В критическом слое построены асимптотические разложения для ¥п1 и кп1. Комплексное волновое число определяется

Учитывая формулу (7), получаем, что короткие спиральные волны движутся с характерной скоростью стационарного потока и.

Функция Нп т в (7) локализована в области пограничного слоя вблизи оболочки

и исчезает при выходе из этого слоя. Разлагаем функцию Нп т в асимптотический ряд.

Главный член Н этого разложения определяется из неоднородной краевой задачи:

Н = 0 (¿1 = 0), Н = 0 (¿1 = да).

Оценивая правую часть в полученном уравнении, находим порядок невязки для т = 1, п > 1:

отметим, что невязка убывает с уменьшением параметра , с ростом индекса п и при выходе из области пограничного слоя.

В окрестности оболочки при г = 1 + О(1) для конечных значений функции Иг 0 и при больших значениях параметра Я приведем формулу для оценки азимутальной

формулой

у

р

компоненты скорости П

(я ^ да). Итак, короткие спиральные волны

сосредоточены вне пристеночного пограничного слоя.

Зависимость окружной компоненты скорости п01 от времени рассчитывалась по формуле (7) для первой моды т = п = 1 при фиксированных значениях параметров и пространственных координат. Правая часть формулы (7) умножалась на постоянную с = сг + гс{. Расчеты показали, что в зависимости от значений коэффициентов сг, с;

возможны различные варианты закрученных потоков крови в систолу. Например, может возникать однонаправленное вращение. Однако это направление вращения может изменяться на противоположное в течение небольшого промежутка времени. Этот промежуток может находиться вблизи любой временной точки в период систолы, а его длительность мала по сравнению с ее продолжительностью. Численные расчеты подтверждают опубликованные в [7] результаты экспериментальных исследований, проведенных на восходящем отделе аорты собаки.

Замечание

При отсутствии стационарного потока У0 = 0 декремент затухания спиральных

волн имеет большое значение порядка О

и спиральные волны затухают сразу

вблизи входа в аорту. В случае У0 ^0 декремент затухания мал и имеет порядок О(8У) .

Итак, механизмом переноса коротких спиральных волн является стационарный поток, так как при У0 = 0 перенос спиральных волн полем длинных продольных волн

незначителен. Для аорты собаки спиральные волны в течение систолы могут переноситься на расстояние от двух до семи сантиметров в зависимости от величины скорости стационарного потока, изменяющейся в пределах от 10 до 40 см/с [8]. Систола - одна из фаз сердечного цикла - сокращение. Во время систолы кровь нагнетается в артериальную систему. В расчетах полагалось, что систола составляет одну треть сердечного цикла для аорты собаки [2].

Квазистационарные моды

Спиральные моды (7) изменяют направление вращения жидкости на обратное либо по сечению, либо со временем в течение сердечного цикла. Построим асимптотические разложения квазистационарных мод, которые в главном приближении не зависят от времени, причем зависимость от времени проявляется только в высших приближениях. Малые параметры 8к и 8у свяжем по формуле 8к = , где с8 = О(1).

Для аорты собаки [8] приведем значение с8« 3,44 . Вне области пограничного слоя асимптотические разложения квазистационарных мод представим в виде

= 8^0 + 82Щ + ••• , (8к ^0). (8)

Подставляем разложение (8) в уравнение (2) и учитываем асимптотические ряды (3) для длинноволновых компонент вектора скорости , уг . В результате находим, что д ^

---0 = 0. Итак, главное приближение ряда (8) не зависит от времени ^0 = w0(г, .

д Ї

Далее выводим уравнение для функции w1 в (8), из которого легко находится общее решение. Выполняя условие периодичности по времени для функции w1, выводим уравнение

+ 1 ап - W2L - ,8КГЛг)^ = 0.

д г г д г г д

2

1

г

0 0,25 0,5 0,75

1

Рис. 2. Квазистационарная мода

Краевое условие на стенке цилиндра для функции п0 получаем, подставляя (8) в динамические уравнения оболочки. В результате выводим формулу п0 = 0 (г = 1). На оси цилиндра выполняется условие п0 = 0 (г = 0). Решение краевой задачи для

функции п0 представим в виде формулы п0 = Ш0(г)ехр . В последней формуле

учтено, что медленная переменная выражена через переменную 2 , и введено число

и а

Рейнольдса Яе = ——. Параметр X находится численно путем решения задачи на

V

собственные значения. В результате находим счетное число собственных чисел X = Хк

и собственных функций = п0к(г) (к = 1,2,3,...). Отметим, что для первой моды

X = 21,382, а для второй X2 = 74,779. На рис. 2 изображена зависимость от

радиальной координаты первой собственной функции п01(г) по норме пространства Ь2

[0,1]. Функции п0к (г) при к > 2 обращаются в нуль (к -1)-раз внутри промежутка

[0, 1] и принимают как положительные, так отрицательные значения. Первая мода главного приближения ряда (8) не изменяет направления вращения ни со временем, ни по сечению и совместно с короткими и длинными спиральными волнами описывает винтовые течения жидкости в аорте в течение сердечного цикла. Учет пограничного слоя приводит к разложению

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

В работе построены асимптотические разложения коротких спиральных волн и квазистационарных спиральных мод, описывающих течения жидкости в аорте, которая моделируется цилиндром, ограниченным тонкой упругой оболочкой. Эти волны, в отличие от длинных спиральных волн [10, 11], заполняют все поперечное сечение цилиндра и слабо зависят от упругих свойств оболочки. Механизмом переноса коротких спиральных волн является стационарный поток, который переносит эти волны на небольшое расстояние в течение сердечного цикла в зависимости от скорости потока.

^ =8 к^( г, 2і) +8 2 ^(г, ^ і ) +8 н 0 + О(8к)

Выводы

Благодарность

Работа выполнена при финансовой поддержке гранта РФФИ № 12-01-00582-а.

Список литературы

1. Багаев С.Н., Захаров В. А., Орлов В. А. О необходимости винтового движения крови // Российский журнал биомеханики. - 2002. - Т. 6, № 4. - С. 30-51.

2. Богаченко С.Е., Устинов Ю.А. Модель движения крови в артериальном сосуде во время систолы и анализ напряженного состояния стенки с учетом винтовой анизотропии // Российский журнал биомеханики. - 2009. - Т. 13, № 1. - С. 29-42.

3. Бокерия Л. А., Городков А.Ю., Кикнадзе Г.И., Николаев Д. А. Анализ поля скоростей закрученного потока крови в аорте на основании 3D-картирования с помощью МР-велосиметрии // Бюллетень НЦССХ им. А.Н. Бакулева Российской академии медицинских наук. - 2003. - Т. 4, № 9. - С. 70-74.

4. Громека И.С. О скорости распространения волнообразного движения жидкостей в упругих трубках // Собрание сочинений. - М.: Изд. АН СССР, 1952. - С. 172-183.

5. Доброва Н.Б., Кузьмина Н.Б., Роева Л.А. Связь анатомических и гидродинамических особенностей сердца в связи с его насосной функцией // Вестник АМН СССР. - 1974. - № 6. - С. 22-31.

6. Кикнадзе Г.И., Краснов Ю.К. Эволюция смерчеобразных течений вязкой жидкости // Доклады АН СССР. - 1986. - Т. 290, № 6. - С. 1315-1319.

7. Орловский П.И., Гриценко В.В., Углов Ф.Г. Следует ли учитывать наличие закрученного потока крови в левом желудочке сердца и аорте при конструировании искусственных клапанов? // Вестник хирургии. - 1998. - Т. 157, № 1. - С. 10-16.

8. Педли Т. Гидродинамика крупных кровеносных сосудов. - М.: Мир, 1983. - 400 с.

9. Пуриня Б.А., Касьянов В.А. Биомеханика крупных кровеносных сосудов. - Рига: Зинатне, 1980. -260 с.

10. Углов Ф.Г., Зубцовский В.Н., Большаков О.П. Рельеф внутренней поверхности стенки левого желудочка сердца в различной фазе сердечной деятельности и его функциональное значение // Вестник хирургии. - 1984. - № 3. - С. 3-9.

11. Устинов Ю.А. Модель винтового пульсового движения крови в артериальных сосудах // ДАН. -2004. - Т. 398, № 3. - С. 71-76.

12. Фатенков В.Н. Биомеханика сердца в эксперименте и клинике. - М.: Медицина, 1990. - 160 с.

13. Marinelli R.A., Penney D.G., Marinelli W.A. Rotary motion in the heart and blood vessels // J. App. Cardiology. - 1991. - Vol. 6, No. 2. - P. 421-431.

14. Zakharov V.N. The new conception of blood circulation mechanics // Cardiovascular Engineering. - 1998. -Vol. 3, No. 2. - P. 100-104.

INTRA-AORTIC SPIRAL WAVES MODELLING

V.A. Batischev, D.S. Petrovskaya, Y.A. Ustinov (Rostov-on-Don, Russia)

The research presents two types of spiral waves in the aorta - the short spiral waves and the quasi-stationary modes. The study was conducted on the basis of the Navier-Stokes equations and dynamic equations of a thin elastic membrane. The aorta is simulated by a circular cylinder bounded by the isotropic membrane. It is demonstrated that short spiral waves fill in the entire cross-section of the cylinder, as opposed to the long spiral waves, which are located in the boundary layer near the membrane. The elastic properties of the membrane have little effect on the short waves. Part of the spiral modes is located in the critical layer near the axis of the cylinder. The steady-state flux is the transfer mechanism both for the short spiral waves and quasi-stationary modes.

Key words: spiral waves, flow, aorta, elastic membrane, viscosity, boundary layer, quasi-stationary modes.

Получено 18 февраля 2012