Ф. Х. Тазюков, Х. А. Халаф, К. М. Алиев,
Р. С. Шайхетдинова
МОДЕЛИРОВАНИЕ ТЕЧЕНИЙ НЕНЬЮТОНОВСКИХ ЖИДКОСТЕЙ В КАНАЛАХ, СНАБЖЕННЫХ ЗАПОРНЫМ КЛАПАНОМ
Ключевые слова: течение неньютоновской жидкости, канал с запорной арматурой, численное моделирование.
В настоящей работе представлено численное моделирование ламинарного течения неньютоновской жидкости в плоском канале, включающего в себя запорную конструкцию, состоящую из опорного кольца и подвижного клапана. Исследуемые течения характерны для течений физиологических жидкостей в сосудах, снабженных искусственным клапаном, а также течением различных эмульсий в каналах, содержащих запорную арматуру. В качестве реологического конститутивного соотношения используется модель обобщенной ньютоновской жидкости, предсказывающей эффект аномалии вязкости. Дифференциальные уравнения движения совместно с условиями на границах области течения интегрируются методом контрольного объема с использованием алгоритма SIMPLE. Показаны области течения, содержащие пиковые значения напряжений, приводящие к разрушению структуры жидкостей.
Key words: Non-Newtonian liquid flows, channel with disk-type prosthetic heart valve, numerical
simulation.
The present work presents the numerical simulation of the laminar flow through a duct of rectangular section, containing complicated structure. The differential equations that describe the flow were integrated by the Finite Volumes Method, in two dimensions, using SIMPLE algorithm. The mesh is structured, with rectangular volumes. Flows of physiological liquids through a disc-type prosthetic heart valve is investigated. Generalized Newtonian fluid flow is simulated using control volume finite difference method. Peacs of stresses leadind to failure of liquid structure is showed.
Математическая постановка задачи
В последнее время, для определения сдвиговых и нормальных напряжений в жидкости вблизи запорного устройства и улучшении конструкции запорного устройства, математическое моделирование становится все более популярным, особенно для физиологических течений. Застойные области течений, связанные с циркуляционными течениями вблизи клапанов и захватом различных частиц способно приводить к различным местным сужениям. Вблизи стенок канала и на границе между основным потоком и застойной областью могут возникать большие градиенты напряжений, приводящие к разрешению структуры жидкости. Проблемы образования застойных течений, и как следствие образование тромбов, вблизи искусственных сердечных клапанов являются следствием конструктив-
ных недостатков запорных устройств. Вследствии сказанного, моделирование течений в каналах, содержащих запорную арматуру, становится весьма актуальным [1], [3], [4].
В представленной работе исследуется стационарное течение обобщенной ньютоновской жидкости в плоском канале, содержащем полностью открытый клапан дискового типа для различных значений числа Рейнольдса (рис.1). Условиями на границах области течения являются условия прилипания на твердых стенках. Сами стенки являются гладкими, недеформируемыми и непроницаемыми. Во входном сечении задается скорость втекания. В выходном сечении задается установившийся режим течения. На оси канала задаются условия симметрии. Относительные размеры запорного устройства приведены на рис.2.
Уравнения движения обобщенной ньютоновской жидкости записывается в виде
ари^ дру)^0
дх ду
д рии ) д руи) дх су
д риу ) д £>\Л/ ) дх ду
до д ди д ди
— + — и— +— и — дх дх дх ду ду
до д дм д дм
----1----------и— н----и—
ду дх дх ду ду
(1)
(2)
(3)
Рис. 1 - Схема течения в канале
Рис. 2 - Относительные размеры дискового запорного устройства
Для обобщенной ньютоновской жидкости [6]:
ди дм ди дм
ту = 2и—, т. = 2и— , ту..=и— + —
^ н дх УУ н ду ху н ду дх
(4)
где
М = к
М/2
1 Зи 2 Зу 2 Зи Зу 2
где 2 + д )+ ^+ д ^ - вт°р°й инвариант тензора скоростей деформа-
ций; К - коэффициент консистенции; П - показатель неньютоновости; и и V компоненты скорости в х и у направлениях соответственно; р - плотность жидкости, р - давление, р - эффективная вязкость.
Граничные условия записываются в виде
х = 0 :и = цп,у = 0 х = Ц : ди/ дх = 0 , V = 0
(6)
у = 0 : и = 0 , V = 0 у = Н : Зи / Зу = 0 , V = 0
В безразмерном виде уравнения (1)-(3) запишутся в виде
¿но+сМ/_)=0 дх ду
д0) д0) Зр Зтхх Зтху
Ке <|1^Г+''^)р^дГ+д' + д № (|^и+у^-))=-.5Р+ ^ + ^
дх ду ду дх ду ’
где Ке = р1_п\/2'п / К - обобщенное число Рейнольдса; К 5"1 - характерная вязкость;
- характерные длина и скорость соответственно.
Метод решения
Для численного решения уравнений (1)-(5) с граничными условиями (6) использован метод контрольного объема, включающий в себя интегрирование уравнений движение по контрольному объему.
Для дискретизации области течения используется шахматная сетка. Область течения делится на контрольные объемы, отределяемые штриховой линией (рис. 3). Давление рассчитывается в узловых точках ( • ) определяемых как пересечение сплошных линий и отмеченных заглавными буквами Р, ', Е, N и Б. Компонента скорости и вычисляется на сторонах контрольного объема (КО), помеченных строчными буквами и . Компонен-
та скорости V вычисляется на сторонах контрольного объема (КО), помеченных строчными буквами и . Величины и Ау являются шагом сетки в и направлениях соответственно.
Уравнение сплошности в дискретной форме записывается как
МО - ^и1^у+[^И- = 0. .......
(7)
Рис. 3 - Контрольные объемы для различных переменных
Дискретные формы для уравнения движения в проекциях на оси координат
а +е = I апь+пь + Ь п + рР - РЕ ^у
(8)
и и и д)
Д|. Е W S N ’ е Д|.
ehe =XanV,hnb+bE)+ pN-PP^
где
где
(О
(9)
At
At
где Ре , Р°, и° и У° относятся к известным значениям в момент времени {, все остальные величины являются неизвестными величинами определяемыми в момент времени .
Коэффициенты с индексами ( ) и ( ) относятся к величинам в и направлениях соот-
ветственно. Коэффициенты а
nb 5 anb
относятся к коэффициентам а^ \ aw
о)
у
яЕ , aw , aN и аs ' , которые определяют взаимное конвективно-диффузионное влияние потока на сторонах и и V контрольного объема соответственно. Значения этих коэффициентов определяются на основе степенной схемы [2]. Компоненты скорости unb и vnb являются значениями скоростей в соседних узловых точках.
Численная процедура
Для решения системы дискретных уравнений (7)-(9), используется алгоритм SIMPLE algorithm [2]-[3], являющийся по существу процедурой предиктор-корректор для расчета давления на смещенной сетке. Для инициализации алгоритма SIMPLE задается поле давления р *. Уравнения движения в дискретной форме (8) и (9) , для заданного поля давления решаются для получения компонент скорости и
ье =13^ +ЬЕ)+ ^р-РеК
апЕ+п =13^+13^+ ^-P^. .......
(11)
Определим коррекцию давления р ' как разность между уточненным давлением р и выбранным приближенным значением давления р * . Таким образом можно записать
Р = Р*+Р- .......(12)
Аналогично определяется коррекция скорости и ‘ и V ■ как разность между уточненной ( и ) и приближенной скоростями ( и )
u = u*+u- ... (13)
.
(14)
Вычитание уравнений (10) и (11) из уравнений (8) and (9) приводит к
аЕ+е-и*е ^ ХЭпЬ+пЬ -u *пЬ Ир-Р*Р I ^-Р*еК .......
(15)
a„n)i -v*n ^la^+ta -v*nb -P*n - ^-Р*р^х .......
(16)
Используя формулы для коррекции давления и скоростей (12)-(14), уравнения (15) и (16) могут быть записаны в виде
3^=23^+ рр'-РЕ^у
аеЕ *V= 1ап^пь + pN-PP^x
(17)
(18)
Согласно алгоритму SIMPLE, для упрощения уравнений (17) и (18), слагаемыми и Ха^^пь можно пренебречь. Следовательно, для коррекции скоростей можно получить следующие выражения
и.-=^п рр--РЕ '
(19)
е (j ) ЧР “ 1 Е
vn"_ я v ) PlsT- рр"
Подстановка уравнений (19) и (20) в (13) и (14) приводит к
ue - u *e + „ п ) ^P+ PE
(21)
8 ae^) E
Ax / •
........
(22)
Аналогично
Ay
Uw U w + „W) І+ - PP
(23)
a(J) 'VW P dw
Ax і • SV%+-^i tfp'-Ps'.
(24)
Подстановка выражений (21)-(24) в дискретное уравнение сплошности (7), приводит к написанию уравнения для определения коррекции давления, играющее ключевую роль в алгоритме SIMPLE
арРр" аЕРЕ+ awPw + aNPN+ aSPS+ [ *W - u ^ ^ v ^ - ^v * ^x ............
(25)
где ap - aE+aw + aN + as
___ Ду2 _ __ &2 ___ Ax2 ___ ^2
aE “ Pe a ; 3W “ Pw _ ; 3N “ Pn _ ; 3S “ Ps a
“w n “s
Процедура алгоритма SIMPLE показана на блок-схеме, представленной на рис. 4.
Таким образом, задача сводится к решению системы линейных уравнений вида
арФр=аЕФЕ+аУуЧ%+амФм+азФз+Ь ......
(26)
где
aE=DeA |^е| )*" tt-Fe,0J , aw = DWA |l^w| )*-[[fw,0]]
aN=DnAp|H-Fn,0j , as=DsA^s^[-Fs,^
aP - л + ^ , ^ - ScA)^\y + арФ0 , ap “ aE + aw + aN + as + ap - SpAxAy
Здесь PP и фр являются известными величинами, относящимися к временному слою t.
Рис. 4 - Алгоритм SIMPLE
Потоки Ре , ^ , Рп и Р5 , а также соответствующие коэффициенты диффузии определяются как
Р^^Ду , РИ = ¿и^Ду , Рп = ^Дх (^Дх ,
С>е = ГеДу/Х , 0И=ГИДу/Х , О п = ГпДх/у , О3=Г5Дх/
Числа Пекле
Pe-Fe/De , Pw-^/Dw , ^-^/Dn , Ре -Fe/De
Функция А|*|Н°, Х-ОЛИПЛ определяет степенную схему аппроксимации задачи.
Результаты моделирования
На рис. 5 представлено распределение функции тока для случая модельного дискового клапана при заданном значении числа Рейнольдса ЇЧе=90 и различных значений показателя неньютоновости при одинаковом расстоянии от кольца до запорного диска. На рисунках хорошо видна циркуляционная область, образующаяся вниз по потоку за выступом кольца и за диском.
1
>0.5
Re=90, n=0.8, Stream Function
1
>0.5
>0.5
1 2 3 4 5 6 X Re=90, n=1.0, Pressure P
■H n
1 2 3 4 5 6 X Re=90, n=1.2, Stream Function
їаШШШЖ
1 2
3
X
4 5
1
6
Рис. 5 - Изолинии функции тока при различных значениях показателя неньютоново-сти
Большинство физиологических жидкостей обладает свойством псевдопластичности, то есть уменьшения вязкости с ростом сдвиговых напряжений. Можно видеть, что для псевдопластиков, область циркуляционных течений значительно больше. Следовательно, более густая кровь обладающая меньшим значением показателя неньютоновости более активно способствует образованию различных бляшек (атером) на стенках сосудов. Тромбы, образующиеся вокруг диска, могут отрываться и приводить к повреждению запорного устройства клапана.
Известно [1], что касательные напряжения на стенках сосудов также влияют на образование атером. Имеются также сведения [4], что возникновение сформированных атеросклеротических бляшек в основном происходит вблизи стенок подвергающимся повышенным или пониженным сдвиговым напряжениям. На рис.6 приведены данные по распределению пристенного сдвигового напряжения для различных значений показателя неньютоновости. Видно, что вблизи запорного диска на стенках сосудов образуются пики
сдвиговых напряжений, способные привести к разрушению структуры жидкости и образованию атеросклеротических бляшек.
Рис. 6 - Пристенное сдвиговое напряжение при течении обобщенной ньютоновской жидкости при Рв=90
Из рисунка следует, что при течении псевдопластиков (n<1) на стенках сосудов пра-кически не образуются пики напряжений, тогда как для дилатантных жидкостей (n>1) пики напряжений становятся значительными.
Выводы
В данной работе проведены исследования влияния аномалии вязкости на размеры отрыва потока (циркуляционной области) в канале, снабженным искусственным дисковым запорным устройством, а также на величину пристенного сдвигового напряжения вблизи данного типа запорного устройства. Получено, что эффект аномалии вязкости существенно влияет на все характеристики течения, в том числе как на размеры области циркуляции, так и на образование пиков пристенных сдвиговых напряжений. Проведение численных экспериментов может позволить лучше разобраться в работе запорной арматуры и уменьшить негативные эффекты, возникающие при течении физиологических жидкостей непосредственно связанные с размерами циркуляционных течений и величинами пристенного сдвигового напряжения.
Литература
1. Peskin, C.S. The fluid dynamics of heart valves: experimental, theoretical and computational methods. Ann. Rev. fluid. Mech. V.14, 1982, p.235-259.
2. Patankar, S. V. Numerical Heat Transfer and Fluid Flow, Hemispher, New York, 1980.
3. Х.А. Халаф, К.М. Алиев, Ф.Х. Тазюков Ламинарное течение обобщенной ньютоновской жидкости в канале содержащем открытый клапан. VII школа-семинар молодых ученых и специалистов «Проблемы теплообмена и гидродинамики в энергомашиностроении», Казань, 2010, 238-241 с
4. Thalassoudis, K. and Mazumdar, J, Mathematical Model fo Turbulent Blood Flow through a Disk-type Prosthetic Heart Valve, Medical and Biological Engineering and Computing, Vol. 22, Pp.529-536, 1984.
© Ф. Х. Тазюков - д-р техн. наук, ст. науч. сотр. каф. технологии конструкционных материалов КГТУ, [email protected]; Х. А. Халаф - асп. той же кафедры; К. М. Алиев - лаб. той же кафедры; Р. С. Шайхет-динова- ст. препод. каф. технологии конструкционных материалов КГТУ, [email protected].