CC BY
5
ISSN 1026-2237 BULLETIN OF HIGHER EDUCATIONAL INSTITUTIONS. NORTH CA UCASUS REGION. NATURAL SCIENCE. 2022. No. 4-1
Научная статья УДК 532.5, 556
doi: 10.18522/1026-2237-2022-4-1-42-54
НАЧАЛЬНО-КРАЕВАЯ ЗАДАЧА О ТЕЧЕНИИ В КРОВЕНОСНОМ СОСУДЕ Наталья Михайловна Полякова
Южный федеральный университет, Ростов-на-Дону, Россия nmzhukova@sfedu.ru
Аннотация. Построено решение начально-краевой задачи для системы двух квазилинейных гиперболических уравнений, описывающей вращательно симметричное безвихревое течение идеальной несжимаемой жидкости в бесконечной цилиндрической области (труба, кровеносный сосуд). Считается, что боковая поверхность области - свободная (мягкая) граница, на которой выполнено кинематическое условие. Динамическое условие на границе моделируется некоторым конститутивным соотношением - зависимостью давления P от площади сечения S цилиндрической области. Для исследования течения использована асимптотическая модель, построенная в работе на основе теории мелкой воды. При конструировании решения начально-краевой, а также начальной и краевой задач для системы двух квазилинейных уравнений в частных производных первого порядка использован метод годографа на основе закона сохранения. Для исследования выбран вариант полиномиального конститутивного соотношения: P~S2e (в>0). В случае начальной задачи с данными в начальный момент времени указана функция Римана - Грина, позволяющая сконструировать неявное аналитическое решение, и приведен численный алгоритм построения явного решения исходной задачи. Основное внимание фокусируется на частном случае P~S2, для которого все необходимые для построения решения формулы записываются в явном виде. При P~S2 для начально-краевой задачи рассмотрено несколько вариантов начальных и краевых условий, которые позволяют детально проследить эволюцию решения. Приведены численные расчеты, демонстрирующие движение ударных волн и волновых фронтов. Полученные результаты на практике можно использовать для описания течения в кровеносных сосудах, а также как надежные тесты для проверки вычислительных алгоритмов, предназначенных для решения подобных задач, в частности систем квазилинейных гиперболических уравнений.
Ключевые слова: квазилинейные гиперболические уравнения, начально-краевая задача, сильные и слабые разрывы, течение крови
Для цитирования: Полякова Н.М. Начально-краевая задача о течении в кровеносном сосуде // Известия вузов. Северо-Кавказский регион. Естественные науки. 2022. № 4-1. С. 42-54.
Статья опубликована на условиях лицензии Creative Commons Attribution 4.0 International (CC-BY4.0).
Original article
INITIAL BOUNDARY VALUE PROBLEM OF THE FLOW IN A BLOOD VESSEL Natalia M. Polyakova
Southern Federal University, Rostov-on-Don, Russia nmzhukova@sfedu.ru
Abstract. The solution of the initial boundary value problem for a system of two quasi-linear hyperbolic equations, which describes a rotationally symmetric vortex-free flow of a viscous incompressible fluid in an infinite cylindrical domain (pipe, blood vessel) is constructed. It is assumed that the side surface of the domain is a free (soft, compliant) boundary on which the kinematic condition is set. The dynamic condition at the boundary is
© Полякова Н.М., 2022
modeled by some constitutive relation - the dependence of the pressure P on the cross-sectional area S of the cylindrical domain. To study the flow an asymptotic model based on the shallow water theory is used. When constructing the solution of the initial-boundary (as well as initial and boundary) problem for a system of two quasilinear partial differential equations of the first order, the hodograph method based on the conservation law is used. A variant of the polynomial constitutive relation is chosen for the study: P~S2e (в>0). In the case of the initial data problem which specified at the initial time the Riemann-Green function, allowing to construct an implicit analytical solution, and an algorithm (numerical) for constructing an explicit solution of the original problem are given. The main attention is focused on the special case P~S2, for which all the formulas necessary to construct a solution are written in explicit form. For P~S2, several variants of conditions (initial and boundary value) are considered for the initial boundary value problem, which allow us to trace the evolution of the solution in detail. Numerical calculations demonstrating the motion of shock waves and wave fronts are given. The results obtained, in practice, can be used to describe flows in blood vessels, as well as reliable tests for testing computational algorithms intended to solve such problems, in particular, systems of quasi-linear hyperbolic equations.
Keywords: quasi-linear hyperbolic equations, initial boundary value problem, strong and weak discontinuities, blood flow
For citation: Polyakova N.M. Initial Boundary Value Problem of the Flow in a Blood Vessel. Bulletin of Higher Educational Institutions. North Caucasus Region. Natural Science. 2022;(4-1):42-54. (In Russ.).
This is an open access article distributed under the terms of Creative Commons Attribution 4.0 International License (CC-BY 4.0).
Введение
Моделированию течения крови в крупных кровеносных сосудах (артерии) посвящено значительное количество работ [1-5]. Теоретические, экспериментальные и клинические принципы течения крови отражены, в частности, в монографиях [6, 7]. Одна из простейших и популярных моделей представляет собой систему двух квазилинейных уравнений гиперболического типа. По-видимому, впервые такая модель, построенная при помощи метода осреднения, предложена в [1] и затем существенно развита и исследована многими авторами. Заметим, что практически такая же модель построена в [6, с. 72, 73] на основе феноменологических соображений. Наиболее полно способ построения модели методом осреднения представлен в работе [2], в которой, помимо прочего, проведен анализ гиперболических уравнений, даны численные схемы их решений и имеются результаты расчетов. Следует сказать, что в работах [1-5] одна из самых важных проблем методов осреднения - замыкание уравнений - решается крайне спорно. При выводе уравнения движения использован коррекционный коэффициент, связывающий квадрат средней скорости со средним квадрата скорости, необоснованно считающийся постоянным.
В предлагаемой работе применен иной, отличный от метода осреднения, способ конструирования уравнений - лагранжев подход теории мелкой воды. Для безвихревого течения жидкости построена асимптотическая модель. При помощи разложения по малому параметру (относительный размер сосуда) получено точное решение для функции тока с последующим проектированием решения на боковую поверхность области.
Другой важной проблемой моделирования, возникающей как для модели осреднения, так и для предлагаемой модели, является выбор конститутивного соотношения, связывающего давление Р и площадь S сечения сосуда, - динамического условия на границе области. Как правило, зависимость Р = P(S) выбирается на основе полуэмпирических гипотез [1-6]. Наиболее часто используется полиномиальная зависимость: P~S2P (р > 0).
Очевидно, что с точки зрения гидродинамики вращательно симметричное течение в кровеносном сосуде представляет собой течение в некоторой цилиндрической области (труба, струя). Боковую границу цилиндрической области (стенка трубы) можно считать, например, свободной (мягкая, податливая) и задавать на границе кинематические и динамические условия. Именно такому случаю и соответствует задание Р = P(S) в качестве динамического условия. Можно также учитывать поверхностное натяжение на свободной границе, что более приемлемо при рассмотрении струйных течений. Более реалистично рассматривать упругие или вязко-упругие стенки трубы. Такие варианты моделей [6, 8-10] являются существенно более сложными (вместо алгебраической зависимости Р = P(S) возникают некоторые дифференциальные уравнения), их исследование возможно только численными методами и в представленной работе
не рассматриваются. В связи с другими моделями упомянем [11], в которой приведен обзор моделей, а также [12], где рассмотрена модель, учитывающая наличие источников. Заметим, что во всех моделях течение жидкости предполагается сравнительно медленным. Задача о стационарном турбулентном течении в гофрированной трубе, стенки которой имеют податливые (мягкие) участки, рассмотрена, например, в [13].
С математической точки зрения представленные результаты являются описанием метода построения решения начальной, краевой и начально-краевой задач для системы двух уравнений гиперболического типа. Предлагаемый метод решения (вариант метода годографа на основе закона сохранения [14, 15]) представляет собой численно-аналитический метод, который выгодно отличается от прямых численных методов, применяемых, например, в [1-8] непосредственно к решению задачи о течении крови по сосудам. Следует сказать, что в [1-8] большое внимание уделяется теоретическому качественному анализу свойств гиперболических уравнений, хотя полное решение, аналогичное приведенным, например, в [14] для иных задач, описываемых гиперболическими уравнениями, отсутствует.
В представленной работе выбрана модель с конститутивным соотношением P~S2, что позволило существенно уменьшить громоздкость приводимых результатов, не искажая общей схемы построения решения и описания метода. Более того, указанный выбор позволил расщепить систему двух квазилинейных уравнений на отдельные, хорошо известные и исследованные уравнения Хопфа, которые связаны между собой лишь начальными данными. Это, в частности, дало возможность сравнить результаты, имеющиеся для уравнения Хопфа, с представленными в работе решениями и убедиться в эффективности предлагаемого метода.
Постановка задачи
Для моделирования вращательно симметричного течения идеальной несжимаемой жидкости в бесконечной цилиндрической области С используем уравнения Эйлера, записанные в безразмерных переменных в цилиндрической системе координат (г,в,г), полагая отсутствующей азимутальную компоненту скорости
(ru)r + (rw)z = 0, ut + uur + wuz = —pr, wt + uwr + wwz = —pz, (1)
С = {(r,z):0 <r < R(z, t),-<x, < z < +œ). Здесь u,v,w - радиальная, азимутальная (v = 0), осевая компоненты скорости течения v = (u,v,w); р - давление; С - область, в которой происходит течение; г = R(z,t) - уравнение боковой поверхности области С.
Боковую поверхность считаем свободной границей жидкости, на которой выполнено кинематическое
Rt + wRz = и, г = R(z, t) (2)
и динамическое условия. Последнее при отсутствии вязкости и поверхностного натяжения соответствует равенству на границе давления жидкости р и давления окружающей внешней среды Р
р = Р, r = R(z,t). (3)
На оси г = 0 цилиндрической вращательно симметричной области С выполнены естественные условия
и = 0, wr = 0, г = 0. (4)
Предположение о безвихревом характере течения влечет требование обращения в нуль вихря скорости
ш = rotv, rotrv = 0, rotgv = uz — wr = 0, rotzv = 0. (5)
Для полной постановки задачи следует добавить начальные условия, которые будут поставлены далее при решении конкретных задач.
Основная цель работы - построение асимптотического варианта задачи (1)-(5) с учетом различия характерных размеров области С и нахождение решения для различных начальных, краевых и начально-краевых условий.
Заметим, что использование для описания течения вязкой жидкости уравнений Эйлера (1) взамен уравнений Навье - Стокса связано лишь с целью упрощения изложения. При построении асимптотической модели показано, что члены уравнений Навье - Стокса, связанные с учетом вязкости, исчезают.
Асимптотическая модель
При конструировании асимптотической модели используем лагранжев подход, применяемый для построения уравнений типа мелкой воды. Считая, что течение происходит в цилиндрической области, характерный радиус которой много меньше характерного размера в осевом направлении, определяем малый параметр £, характеризующий соотношение размеров области. Решение исходной задачи (1)-(5) разыскиваем в виде рядов по параметру £ с последующим проектированием уравнений на границу области [16, с. 32-39]. Ключевую роль при построении асимптотической модели играет уравнение (5). Требование безвихревого характера течения после замен г ^ £ г, и ^ £и, (£ ^ 0) и введения функции тока ru = —ipz, rw = ~фг с учетом соотношений (4) приводит к задаче
+ *\~-0. (6)
Поиск решения недоопределенной (не хватает краевого условия) задачи (6), которое возможно лишь с точностью до некоторой произвольной функции, в виде формального ряда по степеням £2 приводит к соотношению ^(r,z,t) = (1/2)r2F(z,t) + 0(£2), где F(z,t) -произвольная функция. Радиальная и осевая скорости, вычисленные при помощи функции тока ^ с точностью до членов порядка 0(£2), имеют вид
u(r,z,t) = -1rFz(z,t) + 0(£2),w(r,z,t) = F(z,t) + 0(£2). (7)
В принципе, соотношения (7) строятся в виде бесконечных рядов и могут быть записаны с любой точностью. При этом такие ряды содержат только одну неизвестную функцию F(z, t) и ее всевозможные производные.
Для определения в соотношениях (7) функции F(z,t) используем краевые условия (2), (3) на свободной границе г = R(z,t). Дифференцируя (3) по z, получим Pz — pr(R, z, t)Rz — pz(R, z, t) = 0. Заменяя pr(R, z, t)Rz, pz(R, z, t) при помощи уравнений (1), имеем
Pz + £2Rz™ + ™ = 0. (8)
z z Dt Dt v '
Здесь U(z,t) = u(R,z,t), W(z,t) = w(R,z,t), т. е. компоненты скорости на границе г = R^z^)^ учетом (7) записываются в виде
U(z,t) = —lRFz + 0(£2), W(z,t) = F + 0(£2). (9)
DU/Dt = Ut + WUz, DW/Dt = Wt + WWz - материальные производные на границе, определяемые с учетом условия (2) [17, с. 148-151].
Используя соотношения (9) для подстановки в кинематическое (2) и динамическое (8) условия и сохраняя лишь члены порядка 0(1), получим асимптотическую модель типа классической мелкой воды [16, с. 36], которую записываем в консервативной форме
St + (SW)Z = 0, S = R2, (10)
Wt + (1W2) +Pz = 0, W = F. (11)
Здесь S - площадь сечения цилиндрической области с точностью до множителя п, функция W при £ ^ 0 совпадает с функцией F с точностью до членов порядка 0(£2).
Уравнения (10), (11) дополняем конститутивным соотношением - зависимостью давления внешней среды от радиуса цилиндрической области
P = P(S). (12)
Система квазилинейных уравнений (10)-(12) формально совпадает со многими моделями течения крови, полученными на основе метода осреднения [1, 2, 5]. Однако модель на основе теории мелкой воды является более гибкой в том смысле, что позволяет определить скорость течения жидкости внутри кровеносного сосуда при помощи (7) и значения F( , ) = W( , ) на границе, что, конечно же, невозможно при использовании модели на основе метода осреднения.
Построение решения
Система квазилинейных уравнений (10), (11) для различных конститутивных соотношений P = P (S) достаточно широко исследована численными и аналитическими (качественно) методами, в частности в [1-5], где, помимо построения решения, обсуждаются также различные варианты зависимости P = P(S). Сравнительно часто используется полиномиальная зависимость
S2p
P(S)=1J' P>0' (13)
где fi - некоторая константа (масштабный множитель для давления задан равным единице). Наиболее часто считают, что fi = 1 /4.
Для построения решения задачи (10), (11), (13) используем метод годографа, основанный на законе сохранения, впервые предложенный в [15], а затем модифицированный и развитый в [14]. Систему уравнений (10), (11) записываем при помощи инвариантов Римана R1(z, t), R2(z , t)
Rj + A1Rl = 0, R? + A2R2 = 0, (14)
где A1(R1, R2), A2(R1, R2) - характеристические направления.
Уравнения (14) дополняем условиями
R\z,t)lr = Rl(T), R2(z,t)lr = R2(z), (15)
Г = {(z, t):z = z(t), t = t(r)], (16)
где Г - кусочно-гладкая ориентированная линия, параметризованная при помощи параметра т; Ro(t), Ro(t) - заданные на Г функции.
В случае (13) связь инвариантов Римана с исходными переменными и функции A1(R1,R2), A2(R1,R2) определяются соотношениями [2, 4]
s? s?
R1 = W -j<R2 = W + j, (17)
W = 1(R1 + R2), =^(R2 -R1), (18)
A1,2 = 1±ËR1+1±lR2. (19)
На самом деле и в общем случае Р = P(S) можно указать связь инвариантов Римана с исходными переменными W, S, подобную (17). Однако крайне затруднительно определить явный вид обратных зависимостей вида (18) и связи (19) (явные выражения для A1,2(R1,R2)), играющих важную роль в методе годографа.
Для построения неявного двухпараметрического решения задачи (14)-(16) следует рассмотреть закон сохранения: (рt(R1,R2) + xpz(R1,R2) = 0, выполнить дифференцирование, с учетом (14) получить уравнения для функций (р, гр: ipRi = A1^ri, ipR2 = A2(pR2 и, используя условие разрешимости уравнений, получить задачи для определения (р, гр. Опуская подробности (детальные построения см. в [14]), укажем, что неявное решение, зависящее от двух параметров: а, Ь, а также функции t = t(a, b), z = z(a, b), необходимые для сопоставления переменных (z, t) параметрам (a, b), записываются в форме
R1(z, t) = R!(b), R2(z, t) = Rg(b), 2t(a, b) = t(a) + t(b) - L (^tdt - ^tdz), (20)
' ab
2z(a, b) = z(a) + z(b) - L (ipzdt - (pzdz), (21)
' ab
где z(a), z(b), t(a), t(b) - координаты и моменты времени на кривой Г, соответствующие значениям параметров а, Ь; ГаЬ - часть линии Г с направлением ориентации от параметра а до Ь.
Функция (р1 является решением линейной задачи
Vr1r2 + At^Ri + Bc<Pr2 = 0, (22)
(tf - A1^t)lRi=ri = 1, r1 = R^(b), (^-А2^)1к2=г2 = -1, r2 = R2(a),
12
а' = т-т2, = <23)
Для определения функции гр1 достаточно проинтегрировать одно из уравнений: ipRi = A1^ri , ^r2 =Л2Фк2 .
Аналогичным образом находятся функции ф2, . Задача для их определения формально записывается путем использования в соотношениях (22), (23) взаимозамен t ^ z, (+(р) ^ (±^), 712 1/А1-2.
В [14] показано, что решение задачи (22) тесно связано с построением явного выражения для функции Римана - Грина линейного уравнения в частных производных второго порядка. В случае (13), (19), вычисляя А£, В£ при помощи (23) и решая задачу (22) [14, с. 24-26, 32, 33], опуская громоздкие выкладки, получаем
tfR1 n2| 1 2л _ 2 2F1(a,1-a;1;x) _ (R1-r1)(R2-r2) _ 1+p , .
Ф (K ,K 1Г ,Г ) = |3(r2-r1)a(R2-R1)1-a, X = (R1-R2)(r1-r2), a=2p. (24)
Здесь 2F1 - гипергеометрическая функция Гаусса; для функции ^ указаны как основные аргументы R1, R2, так и параметры г1, г2.
В случае, когда в соотношении (20) или (21) для функций t(a,b), z(a,b) конкретизированы подынтегральные выражения, т. е. в явной форме решена задача (22), имеется способ построения явного решения исходной задачи, по крайней мере численный.
Для определенности считаем, что функция ( a, b) построена в явной форме. Тогда на линиях уровня функции ( a, b) строим явное решение. Пусть имеется параметризованная линия уровня
t( а(т),Ь(т)) = t*, (25)
где т - параметр; t* - момент времени, идентифицирующий линию уровня.
Дифференцируя (25) по параметру , получаем задачу Коши для обыкновенных дифференциальных уравнений (ОДУ)
ar = —tb(a,b), bT = ta(a,b), а\т=0 = a*, b\T=0 = b*, (26)
где (a*, b*) - точка на изолинии (25), соответствующая параметру т = 0, которая может быть найдена путем построения изолиний (25) на плоскости ( a, b) для различных *.
При отсутствии явного выражения (21) координата Z(t) = z(a(r),b(r)) на линии (25) вычисляется путем решения задачи ZT=](a,b), Z\T=0 = Z*, где
J (a, b) = (А2(гг,г2) — A1(r1,r2))tatb, r1 = r1(b),r2 = r2(a).
Окончательно однопараметрическое (параметр ) явное решение задачи на линии уровня (25) записываем в виде
R1(z,t*) = Rl(b(T)), R2(z,t*) = R§(a(r)), t* = t(a*,b*), z = Z(t) = z(a(r),b(r)). (27)
Детально способ конструирования явного решения исходной задачи (14)-(16) и определения начальных значений a*, b*, Z* описан в [14, с. 48-51].
Специальное конститутивное соотношение
С математической точки зрения (и не только) особый интерес представляет частный вариант конститутивного соотношения (13) при ¡3 = 1
P(S) = 1S2. (28)
Дело в том, что в этом случае практически все формулы, необходимые для построения решения, записываются в явной форме. Опуская промежуточные вычисления, приводим функции, которые не зависят от переменных R1, R2
^ = ^ tf = t+rl yz = t+rl ^z = 2rhLm (29)
' Г 2 —_r 1 ^ ' r 2 —r r 2 —r 1 ^ r 2 — r1
Они с помощью (20), (21) позволяют в зависимости от контура ГаЬ определить t(a, b), z(a, b).
Начальная и краевая задачи. Формулы (20), (21) получены интегрированием по замкнутому контуру PQM интегрального аналога (p dt — y dz) = 0 закона сохранения yt + px = 0 с
учетом условий (22). Контур PQM состоит из характеристики QM, на которой сохраняется инвариант R1, характеристики PM, на которой сохраняется инвариант R2, и линии PQ = ГаЬ, на которой определены начальные данные (15). Характеристики разных семейств QM, PM пересекаются в точке M = (t,z) (см. более подробно: [14, с. 20-30; 15]).
Условие пересечения характеристик различных семейств в некоторой точке M = ( , ) является важным (в противном случае решение не существует) и влияет на выбор контура ГаЬ. Для начальной (данные при t = 0) и краевой (данные при z = 0) задач с учетом R2 > R1 (см. (17)) выбор контура очевиден. Приведем соотношения для t(a, b), z(a, b) в случае начальной и краевой задач.
(i) Контур ГаЬ = Г = {z е [a, b],t = 0,b > a} соответствует начальной задаче. Тогда
, Ь—а . br2(a)—ar1(b) , ч
t(d,b^ = r2(a)—ri(b)' Z(a,b)= ^-г-ЧЮ . (30)
(ii) Контур Гаь = Г2 = {t е [a, b],z = 0, a > b} соответствует краевой задаче. Тогда
ar2(a)—br1(b) = (a—b)r1(b)r2(a)
r2(a)—r1(b) , ( , ) r2(a)—r1(b)
Обратим внимание на то, что соотношения (30), (31) можно получить и без использования функций (29) и соотношений (20), (21). Дело в том, что при выполнении (28) система (14)
trn ил _ ar (a)—br (ь) „гп ил _ (a—b)r (b)r (а) ИП
L(a,b) r2(a)—rl(b) , z(a,b) r2(a)—rl(^ . (31)
представляет собой два не связанных между собой уравнения Хопфа RU + R1Rl = 0, R2 + R2RI = 0, (A1 = R1,7? = R2) (связь осуществляется через условия (15), (16)). В частности, для того чтобы получить соотношения (30), достаточно записать решение уравнений Хопфа, используя метод характеристик: z = r1(b)t + b, z = r2(a)t + a, R1 = r1(b) = R^(b), R2 = r2(a) = Ro(a). Определяя t, z, получим (30). Аналогично получаем (31), решая краевую задачу при помощи метода характеристик: z = r1(b)(t - b), z = r2(a)(t - а), R1 = r1(b) = Rfcb), R2 = r2(a) = R2(a).
Начально-краевая задача. Функции ^(т), r^(t) заданы на линии Г = {0 < z < œ, 0 < t < œ}. При выборе контура ГаЬ с Г возможны три варианта.
Во-первых, Г = ГаЬ с {0 < z < œ}, что в точности, включая формулы (30), совпадает с начальной задачей (i). Во-вторых, Г2 = ГаЬ с {0 < t < œ}, что в точности, включая формулы (31), совпадает с краевой задачей (ii). Наконец, специфичный именно для начально-краевой задачи кусочно-гладкий контур Г3, для которого точка а принадлежит оси координат z = 0, а точка b находится на оси t = 0.
(iii) Контур ГаЬ = Г3 = {г£ [0, b],t = 0,Ь > 0} U {t £ [0, a],z = 0, а > 0}.
= b+ar2(a) = r2(a) b+ar1(b) (32)
L(a,u) Г2(а)-г1(ьу Z(a,u) ' (u) r2(a)-ri(b). (32)
Вновь, как и ранее, соотношения (32) можно получить, решая относительно и систему z = r1(b)t + b,z = r2(a)(t - a).
В фиксированные моменты времени t > 0 область z > 0 разбивается на три подобласти (рис. 1, 2) с движущимися границами z1(t), z2(t) между подобластями:
(i) = {Z2( t) <Z<œ}, (ii) = {0<Z< Z1( t)}, (iii) = {Z1( t) <Z< Z 2 ( t)}. В зависимости от начально-краевых данных поведение границ может быть различным. В случае несогласованных начальных и краевых данных на границах между зонами = 1( ), z = z2(t) могут быть как сильные разрывы - ударные волны (рис. 1), так и слабые разрывы -волновой фронт (рис. 2) [18, с. 129-135]. В принципе, возможна ситуация, при которой возникают и автомодельные решения - волны разрежения (не представлено на рисунках).
Ударные волны. Рассматриваем ситуацию, когда в начальный момент времени t = 0 имеется сильный разрыв. Уравнения (10), (11) записаны в консервативной форме и являются естественными законами сохранения массы и импульса, что позволяет получить условия Ренкина -Гюгонио на разрыве (приведены для случая (28)) D[S] = [WS]; D[W] = (1/2)[W2 + S2]; D = dz(t)/dt, где [f] = f(z(t) + 0, t) - f(z(t) - 0,t) - величина разрыва функции f на линии разрыва z = z(t); f(z(t) + 0,t), f(z(t) -0,t) - значения функций справа и слева от линии разрыва; D - скорость движения линии разрыва.
С учетом соотношений (17), (18) (при fi = 1) комбинация исходных условий на разрыве приводит к соотношениям для инвариантов Римана D[ftfc] = (1/2)[(Rk)2]. Таким образом, уравнения (14), условия на разрыве и разрывные начальные данные - это задача о распаде произвольного разрыва, начального или возникшего в процессе эволюции, известная как задача Римана для уравнения Хопфа.
Очевидно, что поведение линии разрыва, например, границы z-y (t) между областями (ii), (iii) определяется задачей Коши для ОДУ, когда в начальный момент = 0 имеется сильный разрыв начальных данных (начальные и краевые условия не согласованы)
^- = D1(t), z1(0) = 0, R2 = const ,D1=!(R1(z1(t) + 0,t) +R1(z1(t)-0,t)). (33) Аналогично записывается задача для определения границы 2( ) между областями (iii), (i) -достаточно произвести взаимозамену индексов 1^2.
На рис. 1 для конкретных начальных данных продемонстрирована ситуация, когда границы между подобластями являются ударными волнами.
Напомним, что для некоторых начальных данных ударные волны могут возникать и как результат опрокидывания профиля решения в процессе эволюции решения (не только в начальный момент времени).
Волновые фронты. Для начально-краевой задачи начальные данные в точке z = 0, t = 0 могут, оставаясь непрерывными, иметь слабые разрывы (разрывы производных). Достаточно типичной является ситуация, когда краевые условия при = 0 зависят от , а начальные условия
при = 0 являются постоянными (или наоборот). Способы построения решения подобных задач, известных как волновой фронт (см., в частности, рис. 2), изложены, например, в [18, с. 129-135].
Ввиду того что для данной задачи система уравнений (10), (11), (28) распадается на два независимых уравнения Хопфа, достаточно рассмотреть начально-краевую задачу
Vt + Wz = 0, V\t>Qz=a = G(t), Gt(0) Ф 0, V\z>0,t=0 = G(0). (34)
Здесь V(z, t) - либо R1(z, t), либо R2(z, t); G(t) - либо R1(t), либо R02(t).
Таким образом, предполагается, что в точке t = 0, z = 0 имеется разрыв производной по t (и по в силу уравнения). Как известно, слабые разрывы движутся по характеристике и не могут ни возникнуть, ни исчезнуть в процессе эволюции. Следуя [18, с. 129-135], ищем решение задачи (34) в виде
V(z, t) = Vo + V1 (z)y + ■■■, y>0, V(z, t) = Vo = G(0) >0, y < 0, y = t — q(z). (35)
Осуществляя подстановку (35) в (34), приравнивая члены при одинаковых степенях y и опуская промежуточные построения, получаем q( ) = / V0 и задачу для транспортного уравнения, описывающего поведение слабого разрыва, dV1/dt — V2/V0 = 0, V1(0) = Gt(0), решение которой легко строится. Опуская промежуточные соотношения, приведем лишь окончательный результат. Разрыв производной определяется соотношением [Vz] = Vz\y=+0 — Vz\y=—0 = (t — t0)—1, t0 = V0/Gz(0) и движется по характеристике z = V0t (y = 0). С течением времени t происходит рост величины разрыва по модулю и стремление к бесконечности при t ^ t0 (градиентная катастрофа). Напомним, что в точке z0 = V0t0, в принципе, возможно образование сильного разрыва.
Результаты расчетов
На рис. 1, 2 представлены результаты вычислений в случаях двух вариантов данных для начально-краевой задачи - сильного и слабого начального разрыва.
1-й вариант. Ударные волны. Рассматриваем полубесконечную цилиндрическую область (труба, кровеносный сосуд) постоянного сечения S0, в которой в начальный момент времени t = 0 (z >0) c постоянной скоростью W0 > 0 течет жидкость (кровь). На границе z = 0 (начало трубы) считается, что сечение трубы такое же, как и при z > 0, т. е. S0, а скорость течения W при t > 0 гармонически пульсирует с амплитудой £ и круговой частотой П.
S\t>0,z=0 = S0, W\t>0z=0 = W0(1 + £COS at), (36)
S\z>0, t=0 = S0, W\z>0,t=0 = W0,
S0 = 0,4, W0 = 1,0, £ = 0,1, П = 1.
В начальный момент времени = 0 имеется сильный разрыв скорости, величина которого £W0. При t = +0 возникают два движущихся сильных разрыва инвариантов R1, R2 (и, естественно, S, W, см. (17), (18) при ¡3 = 1). Поведение движущегося разрыва z1(t) определяется путем решения задачи (33). Для определения закона движения разрыва инварианта R1 численно интегрируется задача Коши
^ = 1(r\0) + r1(a(z1, t))), Z1(0) = 0, (37)
где r1(a) = (W0(1 + £ cos Пa) — S0) - задано при помощи краевого условия (36).
На каждом шаге интегрирования задачи Коши функция a(z1, t) вычисляется как решение алгебраического уравнения для характеристики z1 = r1(a)(t — a), приходящей с оси z = 0 на линию разрыва = 1( ).
Аналогично находится закон движения z2(t) разрыва инварианта R2 - достаточно в соотношениях (37) заменить индекс 1 на 2 и (—S0) на S0.
С малой относительной погрешностью (~ 5 • 10—3) в случае (36) имеем z1(t) ^ D11, z2(t) « D21, D1 = 0,65, D2 = 1,45. Решение задачи, т.е. функции R1(z,t), R2(z,t), S(z,t), W(z, t), в момент t = 10 показано на рис. 1.
В момент t = 10 сильные разрывы находятся в точках z-y « 6,5; z2 « 14,5. Точное решение дают для z1(t), z2(t) значения z1(10) = 6,4699; z2(10) = 14,4916. Решение в области (ii) (при 0 < z1(t)) и области (iii) (при z1(t) < z < z2(t)) получаем, используя (27) и интегрируя задачу Коши (26). В случае области (ii) производные ta, tb вычисляем при помощи соотношений (31)
ISSN 1026-2237 BULLETIN OF HIGHER EDUCATIONAL INSTITUTIONS. NORTH CA UCASUS REGION. NATURAL SCIENCE. 2022. No. 4-1
Рис. 1. Ударные волны в момент t = 10 для случая (36). D± = 0,65; D2 = 1,45; Zy = 6,470; z2 = 14,492; t = 10; W0 = 1;
S0 = 0,4; £ = 0,1 / Fig. 1. Shock waves at time t=10 for case (36). Di=0.65; D2=1.45; zi=6.470; Z2=14.492; t=10; Wo=1; So=0.4; e=0.1
Рис. 2. Волновые фронты в момент t = 5,9 для случая (38). = 0,6; D2 = 1,4; zr = 3,54;
z2 = 8,26; W,, = 1;S0 = 0,4; £ = 0,1 / Fig. 2. Wavefronts at time t = 5.9 for case (38). Dr = 0.6; D2 = 1.4; z± = 3.54; z2 = 8.26; W0 = 1;S0 = 0.4; £ = 0.1
для краевой задачи (при Ь = 10 значения а* = Ь* = 10,0). В случае области (ш) производные 1а, вычисляем при помощи соотношений (32) для начально-краевой задачи (при £ = 10 значения а* = 5,663; Ь* = 0,425). В области (^ (при г > г2({)) решение строится при помощи соотношений (30). В случае (36) оно остается постоянным при > 0.
В процессе эволюции решения Кг(г,1), И2(г,1) (и, естественно, Б(г,1), Ш(г, ¿)) возможно опрокидывание профилей функции, что приводит к образованию ударных волн помимо ударной волны, возникающей ввиду разрыва начальных данных при £ = 0. В частности, для рассматриваемых условий (36) в момент ^ = 10,382 в точке = 3,405 возникает ударная волна. Это хорошо видно на рис. 1 - в момент £ = 10 профиль функций в окрестности z ^ 3,5 уже близок к разрывному. Расчеты показывают (не приведены), что при дальнейшей эволюции решения возникает серия движущихся ударных волн.
2-й вариант. Волновые фронты. На рис. 2 показано непрерывное решение, соответствующее волновым фронтам.
Начальные данные имеют вид Б|£>0,г=0 = ^ (38)
>о,2=о = №о(1 + езтт),
5|г>0,£=0 = ^ ^1г>0Х=0 = Шо,
Б0 = 0,4, Ш0 = 1,0, £ = 0,1, и отличаются от данных (36) соотношением (38).
В начальный момент = 0 функции непрерывны, но имеется разрыв производной по времени для функции Ш, величина которого равна £ ■ йШ0, и соответствующие разрывы производных для инвариантов Римана, имеющие величину £ ■ йШ0 (разрыв для функции 5 отсутствует). Способ построения численного решения практически такой же, как и в случае ударных волн. Отличие заключается в том, что нет необходимости решать задачу (37) для определения поведения ударных волн. Волновые фронты движутся по характеристикам гг( Ь) = = = (Шо - 5о)I, 22(1) = ^ = (Шо + Бо)I, Э = 0,6; Э2 = 1,4. В отличие от случая ударных волн, когда Ь), г2(Ь) определяются численно, соотношения для волновых фронтов являются точными.
На рис. 2 показано решение задачи, т. е. функции R1(z, t), R2(z, t), S(z, t), W(z, t), в момент t = 5,9. Момент времени, при котором разрыв производной по z, например для функции R1(z, t), обращается в бесконечность, равен t0 = 6 ( V0 = W0 — S0, Gz(0) = £fiW0). При решении задачи Коши (26) в случае области (ii) a* = b* = 5,9 (при t = 5,9), а для области (iii) - a* = 3,336, b* = 0 (при t = 5,9). Как и в случае ударных волн, в области (i), т. е. при z > z2(t), решение остается постоянным в течение всего времени эволюции.
Момент времени, при котором разрыв производной по z для функции R2(z, t) обращается в бесконечность, равен t0 = 14 (V0 = W0 + S0, Gz(0) = £ÙW0). Но прежде чем производная функции R2(z,t) достигнет бесконечности в момент t = 14, произойдет опрокидывание профилей функций R1(z,t), S(z,t), W(z,t) (R2(z,t) непрерывна). В момент t = 6 возникнет ударная волна в точке z = 3,6. Естественно, при этом изменится скорость движения ударной волны, которая уже не будет двигаться по характеристике. В момент = 14 аналогичное опрокидывание профиля с образованием ударной волны произойдет и для функции R2(z, t).
Решения на линиях уровня функции x(a, b). На рис. 1, 2 представлены явные решения, когда на линиях уровня t(a,b) = t* (см. (25)) определяются функции R1(z,t*), R2(z,t*), S(z,t*), W( , *) (см. (27)), позволяющие проследить зависимость решения от в фиксированный момент времени *. Аналогично строится решение на линиях уровня функции ( a, b) = *. Взамен системы обыкновенных дифференциальных уравнений (26) следует рассматривать систему: aT = —Xb(a,b), bT = xa(a,b) с соответствующими начальными условиями. Производные Xb(a, b), xa(a, b) вычисляются в зависимости от областей (i) — (iii) при помощи соотношений (30)-(32). Заметим, что производные xa ,хь удовлетворяют соотношениям Xa = À1(r1(a,b),r2(a,b))ta, хь = À2(r1(a,b),r2(a,b))tb [14, с. 47, формула (9.6)]. В частности, это означает, что при решении задач достаточно иметь информацию лишь об одной функции - ( a, b) или x( a, b).
В некотором смысле, особенно при проведении экспериментов, решение на линии уровня * = ( a, b) является более важным, чем на * = ( a, b). Такое решение позволяет получать информацию об изменении характеристик течения ( S( *, ), W( *, )) в зависимости от времени t в фиксированной точке области z = z* (место расположения измерительного прибора).
Опуская подробный ход построения решения, который не отличается от описанного в предыдущих разделах, приведем результаты расчетов для начальных данных (36), соответствующих ударным волнам (рис. 3).
Выбираем * = 3,405, т. е. точку в окрестности образования ударной волны в области (ii) (см. соответствующую область на рис. 1). На интервале времени 0 < < 2 в точке = * (расположение датчика, измерительного прибора) фиксируется область (i), фактически течение с постоянной скоростью (рис. 1, 3). Моменты 2 и 1 определены соотношениями: 2 = */ D1 = 2,352, D2 « 1,448, ^ = z*/D1 = 5,300, D1 « 0,642. В момент t = t2 в точке z* фиксируется разрыв инварианта R2 - граница между областями (i) и (iii) (рис. 1, 3).
На интервале времени 2 < < 1 мимо датчика в точке = * проходит область (ii). В момент = 1 в точке * фиксируется разрыв инварианта R1 - граница между областями (ii) и (iii) (рис. 1, 3). Начиная с момента = 2, в точке * фиксируется решение из области (ii).
Напомним, что в момент td = 10,382 в точке zd = 3,405 возникает опрокидывание профиля и образуется ударная волна (рис. 1, окрестность точки zd = 3,405). Скорость движения ударной волны определяется решением задачи (33) для инварианта Римана R1. Движение ударной волны осуществляется примерно с постоянной скоростью D0 « 0,569 по закону z(t) « zd + D0(t — td). Величина D0 < D-y « 0,642. Это означает, что разрыв, сформировавшийся в процессе эволюции, не догонит разрыва инварианта R1, который имелся в начальный момент времени и движется по закону z-y (t) « D11.
Легко показать, что решение R1(z*, t), R2(z*, t), S(z*, t), W(z*, t) в области (ii), т. е., начиная с t = t-y (рис. 1), является периодическим по времени с периодом Т = 2п/П (П = 1). Действительно, используя для области (ii) соотношения (31), (27) и периодичность функций 1( b), 2( a) по аргументам для начальных данных (36), имеем ( a + Т, b + Т) = ( a, b), t(a + T,b + T) = t(a, b) + Т, R1(z*,t) = R1(b) = R^(b + T) = R1(z*,t + T), R2(z*, t) = R0(a) = R0(a + T) = R2(z*, t + T).
Таким образом, в точке zd = 3,405 в моменты времени t = td + 2пт, m = 0,1,..., возникает серия ударных волн, генерируемых гармоническими периодическими краевыми условиями (36) при z = 0 (рис. 3).
Во избежание недоразумений подчеркнем, что результаты, представленные на рис. 3, соответствуют z„ = 3,405. В этом случае хорошо видны нелинейные искажения инварианта Римана R1(zlf,t), но практически незаметны нелинейные отклонения от периодичности для инварианта R2(z,t,t). С течением времени профиль функции R2(z, t) также опрокинется и возникнут ударные волны. Для того чтобы датчик регистрировал ударную волну, возникающую в процессе решения, его следует размещать в точке z„ > zd, так как в противном случае сформировавшаяся ударная волна не достигнет датчика. Если разместить датчик в точке z ^ 25, то график функции R2(z^,t) (S(z*,t), W(z*,t)) будет иметь такую же пилообразную структуру, как и график функции R1(zt, t).
The case of shock waves
Заключение
Представленные результаты показывают, что течение жидкости (кровь) по трубе (кровеносный сосуд) с мягкими свободными стенками определяется в основном нелинейными эффектами переноса. В выбранной модели вязкие свойства жидкости не оказывают влияния на поведение решения. Понятно, что использование более сложных моделей, например с учетом вязкости, приведет лишь к сглаживанию решений, в частности, сильных разрывов. Сравнение результатов с имеющимися экспериментальными данными [7, с. 300-308] подтверждает сделанное утверждение. Конечно, при сравнении следует помнить, что рассматривалось частное конститутивное соотношение (28), а не более общее соотношение (13). Заметим, что при использовании соотношения (28), учитывая, что уравнения модели расщепляются на два независимых уравнения Хопфа, возможен иной путь построения решения - непосредственное решение начально-краевой задачи для уравнений Хопфа [18, с. 62-66] и сравнение такого решения с представленными в работе результатами. Решение модельной задачи с использованием функции Римана - Грина (24), которое не приведено ввиду громоздкости, показало, что качественное поведение решения практически не отличается от представленного в работе - имеются ударные волны, волновые фронты и т.п.
Отметим, что представленный способ решения начально-краевой задачи является численно-аналитическим. В многочисленных работах о течении крови в артериях [1-8] задача решается либо качественно, либо прямыми численными методами (конечно-разностные, методы конечных объемов), что приводит к погрешности аппроксимации, к погрешности решения больших систем уравнений, возникающих при аппроксимации, к практически неустранимым трудностям вычислений в окрестности сильных и слабых разрывов. Напротив, в представленном в работе способе решения точность ограничивается лишь погрешностью решения задачи Коши для ОДУ, которая существенно проще задачи для уравнений в частных производных. Выбор начальных данных в виде (36) и (38) объясняется их типичностью для исходной задачи, простотой представления и трактовки результатов. Без каких-либо изменений, в частности, основных соотношений
0 5 10 15 20 25
Рис. 3. Зависимость значений функций R1, R2, S, W от времени t при фиксированном z = zt = 3,405; t1 = 5,300, t2 = 2,352 - моменты, в которые регистрируются разрывы решения. Случай ударных волн / Fig. 3. Dependence of the values of the functions R1, R2, S, W on time t at a fixed z = z„ = 3.405. t1 = 5.300, t2 = 2.352 - moments at which discontinuities of the solution are registered.
(29)-(32), построение решения возможно для случая цилиндрической области с переменным радиусом, в сужающемся кровеносном сосуде, иными начальными данными и т.п.
Практическая применимость представленных результатов достаточно очевидна - возможность моделирования течения крови по кровеносным сосудам различного профиля и с различными условиями начального кровотока. В частности, представленные на рис. 3 результаты являются моделированием контроля пульса и, более того, восстановлением картины течения в кровеносном сосуде, по крайней мере качественно, по результатам отдельных неинвазивных измерений.
Список источников
1. Barnard A.C.L., Hunt W.A., Timlake W.P., Varley E.A Theory of Fluid Flow in Compliant Tubes // Biophysical J. 1966. Vol. 6, № 6. P. 717-724. Doi: 10.1016/S0006-3495(66)86690-0.
2. Canic S., Kim E.H. Mathematical analysis of the quasilinear effects in a hypernolic model of blood flow through compliant axi-symmetric vessels // Math. Meth. Appl. Sci. 2003. Vol. 26. P. 1161-1186. Doi: 10.1002/mma.407.
3. Anliker M., Rockwell R.L., Ogden E. Nonlinear Analysis of Flow Pulses and Shock Waves in Arteries. Part I: Derivation and Properties of Mathematical Model // ZAMP. 1971. Vol. 22. P. 217-246. Doi: 10.1007/BF01614000.
4. Formaggia L., Nobile F., Quarteroni A., Veneziani A. Multiscale Modelling of the Circulatory System: a Preliminary Analysis // Computing and Visualization in Science. 1999. Vol. 2, № 2. P. 75-83. Doi: 10.1007/s007910050030.
5. Canic S., Li T. Critical thresholds in a quasilinear hyperbolic model of blood flow // Networks and Heterogeneous Media. 2009. Vol. 4, № 3. P. 527-536. Doi: 10.3934/nhm.2009.4.527.
6. Педли Т. Гидродинамика крупных кровеносных сосудов. М.: Мир, 1983. 400 с.
7. Cardiovascular Fluid Dynamics. Ed. Bergel D. H. L.; N.Y.: Academic Press, 1972. Vol. 1. 365 p.
8. Mikelic A., Guidoboni G., Canic S. Fluid-structure interaction in a pre-stressed tube with thick elastic walls I: the stationary Stokes problem // Networks and Heterogeneous Media. 2007. Vol. 2, № 3. P. 397-423. Doi: 10.3934/nhm.2007.2.397.
9. Волобуев А.Н. Течение жидкости в трубах с эластичными стенками // Успехи физ. наук. 1995. Т. 165, № 2. С. 177-186. Doi: 10.3367/UFNr.0165.199502c.0177.
10. Кудряшов Н.А., Синельщиков Д.И., Чернявский И.Л. Нелинейные эволюционные уравнения для описания возмущений в вязкоэластичной трубке // Нелинейная динамика. 2008. Т. 4, № 1. C. 69-86. Doi: 10.20537/nd0801004.
11. Piccioli F., Bertaglia G., Valiani A., Caleffi V. Modeling blood flow in networks of viscoelastic vessels with the 1-D augmented fluid-structure interaction system // J. of Computational Physics. 2022. Vol. 464. 45 p. Doi: 10.1016/j.jcp.2022.111364.
12. Свиридова Н.В., Власенко В.Д. Моделирование гемодинамических процессов сердечно-сосудистой системы на основе данных периферической артериальной пульсации // Мат. биология и биоинформатика. 2014. Т. 9, № 1. С. 195-205. Doi: 10.17537/2014.9.195.
13. Жуков М.Ю., Полякова Н.М., Ширяева Е.В. Квазистационарное турбулентное течение в цилиндрическом канале с неровными стенками // Изв. вузов. Сев.-Кавк. регион. Естеств. науки. 2020. № 1. C. 4-10. Doi: 10.18522/1026-2237-2020-1-4-10.
14. Жуков М.Ю., Ширяева Е.В., Долгих Т.Ф. Метод годографа для решения гиперболических и эллиптических квазилинейных уравнений. Ростов н/Д.: Изд-во ЮФУ, 2015. 126 c.
15. Senashov S.I., Yakhno A. Conservation laws, hodograph transformation and boundary value problems of plane plasticity // SIGMA. 2012. Vol. 8, № 71. 16 p. Doi: 10.3842/sigma.2012.071.
16. ОвсянниковЛ.В., Макаренко Н.И., НалимовВ.И., Ляпидевский В.Ю., ПлотниковП.И., Стурова И.В., Букреев В.И., Владимиров В.А. Нелинейные проблемы теории поверхностных и внутренних волн. Новосибирск: Наука, 1985. 385 с.
17. ЖуковМ.Ю., ШиряеваЕ.В., ПоляковаН.М. Моделирование испарения капли жидкости. Ростов н/Д.: Изд-во ЮФУ, 2015. 208 с.
18. Уизем Дж. Линейные и нелинейные волны. М.: Мир, 1977. 622 с.
References
1. Barnard A. C. L., Hunt W. A., Timlake W. P., Varley E. A. Theory of Fluid Flow in Compliant Tubes. Biophysical Journal. 1966;6(6):717-724, doi: 10.1016/S0006-3495(66)86690-0.
ISSN 1026-2237 BULLETIN OF HIGHER EDUCATIONAL INSTITUTIONS. NORTH CA UCASUS REGION. NATURAL SCIENCE. 2022. No. 4-1
2. Canic S., Kim E. H. Mathematical analysis of the quasilinear effects in a hypernolic model of blood flow through compliant axi-symmetric vessels. Math. Meth. Appl. Sci. 2003;26:1161-1186, doi: 10.1002/mma.407.
3. Anliker M., Rockwell R. L., Ogden E. Nonlinear Analysis of Flow Pulses and Shock Waves in Arteries. Part I: Derivation and Properties of Mathematical Model. ZAMP. 1971;22:217-246, doi: 10.1007/BF01614000.
4. Formaggia L., Nobile F., Quarteroni A., Veneziani A. Multiscale Modelling of the Circulatory System: a Preliminary Analysis. Computing and Visualization in Science. 1999;2(2):75-83, doi: 10.1007/s007910050030.
5. Canic S., Li T. Critical thresholds in a quasilinear hyperbolic model of blood flow. Networks and Heterogeneous Media. 2009;4(3):527-536, doi: 10.3934/nhm.2009.4.527.
6. Pedli T. The fluid mechanics of large blood vessels. Moscow: Mir Publ.; 1983. 400 p. (In Russ.).
7. Bergel D. H., ed. Cardiovascular Fluid Dynamics. London; New York: Academic Press; 1972;1. 365 p.
8. Mikelic A., Guidoboni G., Canic S. Fluid-structure interaction in a pre-stressed tube with thick elastic walls I: the stationary Stokes problem. Networks and Heterogeneous Media. 2007;2(3):397-423, doi: 10.3934/nhm.2007.2.397.
9. Volobuev A. N. Fluid flow in tubes with elastic walls. Uspekhi fiz. nauk = Physics-Uspekhi. 1995;165(2):177-186, doi: 10.3367/UFNr.0165.199502c.0177. (In Russ.).
10. Kudryashov N. A., Sinel'shchikov D. I., Chernyavskiy I. L. Nonlinear evolution equations for description of perturbations in a viscoelastic tube. Nelineinaya dinamika = Russian Journal of Nonlinear Dynamics. 2008;4(1):69-86, doi: 10.20537/nd0801004. (In Russ.).
11. Piccioli F., Bertaglia G., Valiani A., Caleffi V. Modeling blood flow in networks of viscoelastic vessels with the 1-D augmented fluid-structure interaction system. Journal of Computational Physics. 2022;464:45, doi: 10.1016/j.jcp.2022.111364.
12. Sviridova N. V., Vlasenko V. D. Modelling of hemodynamic processes in cardiovascular system on the base of peripheral arterial pulsation data. Mat. Biologiya i Bioinformatika = Mathematical Biology and Bioin-formatics. 2014;9(1):195-205, doi: 10.17537/2014.9.195. (In Russ.).
13. Zhukov M. Yu., Polyakova N. M., Shiryaeva E.V. Quasi-stationary turbulent flow in cylinder channel with irregular walls. Izv. vuzov. Sev.-Kavk. region. Estestv. nauki = Bulletin of Higher Educational Institutions. North Caucasus Region. Natural Science. 2020;(1):4-10, doi: 10.18522/1026-2237-2020-1-4-10. (In Russ.).
14. Zhukov M. Yu., Shiryaeva E. V., Dolgikh T. F. Hodograph method for solving hyperbolic and elliptic quasilinear equations. Rostov-on-Don: Southern Federal University Press; 2015. 126 p. (In Russ.).
15. Senashov S. I., Yakhno A. Conservation laws, hodograph transformation and boundary value problems of plane plasticity. SIGMA. 2012;8(071). 16 p., doi: 10.3842/sigma.2012.071.
16. Ovsyannikov L.V., Makarenko N.I., Nalimov V.I., Lyapidevsky V.Yu., Plotnikov P.I., Sturova I.V., Bukreev V.I., Vladimirov V.A. Nonlinear problems of theoretical and internal waves. Novosibirsk: Nauka Publ.; 1985. 319 p. (In Russ.).
17. Zhukov M. Yu., Shiryaeva E. V., Polyakova N. M. Modeling the evaporation of a liquid drop. Rostov-on-Don: Southern Federal University Press; 2015. 208 p. (In Russ.).
18. Whitham G.B. Linear and nonlinear waves. Moscow: Mir Publ.; 1977. 622 p. (In Russ.).
Информация об авторе
Н.М. Полякова - ассистент, кафедра вычислительной математики и математической физики, Институт математики, механики и компьютерных наук им. И.И. Воровича.
Information about the author
N.M. Polyakova - Assistant, Department of Computational Mathematics and Mathematical Physics, Vorovich Institute of Mathematics, Mechanics and Computer Science.
Статья поступила в редакцию 01.08.2022; одобрена после рецензирования 14.08.2022; принята к публикации 14.11.2022. The article was submitted 01.08.2022; approved after reviewing 14.08.2022; accepted for publication 14.11.2022.
