CC BY
36
УДК 532.51
Конечно-элементное моделирование неизотермического стационарного течения неньютоновской жидкости в сложных областях
© Ю.И. Димитриенко, Шугуан Ли МГТУ им. Н.Э. Баумана, Москва, 105005, Россия
Метод конечных элементов используется для моделирования неизотермического потока неньютоновских вязких жидкостей в сложных геометриях. Рассмотрена модель Carreau-Yasuda неньютоновской жидкости, в которой зависимость коэффициента вязкости от второго инварианта тензора скоростей деформации имеет степенной вид. Получена вариационная формулировка задачи движения неньютоновской жидкости для плоского случая. Для решения системы уравнений Навье-Стокса применяется итерационный алгоритм Ньютона-Рафсона, а для решения уравнения энергии использован итерационный алгоритм Пикара. Рассмотрена задача о движении полимерной массы в пресс-форме сложного переменного сечения при наличии неравномерного температурного поля. С помощью конечно-элементного моделирования проведен численный анализ влияния различных параметров на движение жидкости и теплопередачу полимерного материала при различных значениях внешнего давления. Показано, что характер движения неньютоновской жидкости существенно зависит от реологических свойств жидкости и характеристик геометрической формы, что необходимо учитывать при технологических процессах переработки пластмасс.
Ключевые слова: метод конечных элементов; неньютоновская жидкость; неизотермическое течение; итерационный алгоритм; сложные геометрические формы, конвекция. .
Введение. Исследованием движений неньютоновских вязких жидкостей занимаются уже почти столетие, в частности их применяют для моделирования движения потоков полимеров в процессе формования изделий из пластмасс [1-3]. Большое значение имеет изучение особенностей движения вязкого полимерного связующего при создании композиционных материалов [4-6], однако ввиду сложности задач моделирования движения неньютоновских жидкостей в пористых системах в настоящее время их изучают лишь в рамках линейно-вязких моделей [7-10]. По сравнению с ньютоновской жидкостью, коэффициент вязкости неньютоновской вязкой жидкости уже не является константой, а зависит, например, от функции инварианта тензора скорости деформации. Кроме того, в процессе литьевого формования полимеров на процесс образования высокомолекулярных продуктов и качество готовой продукции большое влияние оказывает изменение температурного поля, поэтому для моделирования процесса наполнения емко-
стей жидким полимером, как правило, необходимо учитывать температурные эффекты.
Вследствие того, что коэффициент вязкости неньютоновской жидкости зависит от скорости деформаций, численное моделирование процесса ее течения намного сложнее, чем для ньютоновской жидкости [11]. Вследствие зависимости поведения потока неньютоновского вязкого течения от изменения температуры расплава, корректность решения уравнение энергии в значительной степени влияет на итоговый численный результат процесса наполнения емкостей жидким полимером.
В 1965 году Зенкевич и Чунг применили метод конечных элементов (МКЭ) для решения задачи о потенциальном течении жидкости, что стало отправной точкой для решения задачи механики жидкости методом конечных элементов [12,13]. В настоящее время МКЭ наиболее широко применяется для решения различных задач течения неньютоновских жидкостей [14-19]. Различные варианты численных методов решения задач движения неньютоновских жидкостей рассмотрены в [19,20]. Несмотря на совершенствование вычислительной эффективности и точности различных численных методов для неньютоновской жидкости, они, как правило, применимы для относительно простых геометрических форм областей течения, обычно не учитывают тепловых эффектов движения вязкого потока, а также комплексные граничные эффекты. В работах, посвященных методу МКЭ для решения задач течения неньютоновских жидкостей не всегда соблюдается условие ЬББ согласования интерполирующих функций для различных неизвестных -скорости и давления. Поэтому актуальным является дальнейшее развитие численных методов решения задач движения неньютоновских жидкостей, которые
Математическая модель течения полимера в пресс-форме с переменным сечением. Рассматривается течение полимерной жидкости в пресс-форме с переменной площадью поперечного сечения (рис. 1). Обозначим 3 геометрических параметра, характеризующих изменение геометрии сечения: коэффициент сужения 8 = И / Н - где Ь= наименьшая площадь сечения формы, а Н- ее наибольшее значения, и у- угол сужения формы (рис.1).
Общая система уравнений [21] динамики вязкой неньютоновской тепло-проводной жидкости, состоящая из уравнения несжимаемости, уравнений движения и уравнения энергии и записанная в стандартной нотаций индексов [22], имеет следующий вид :
дУ
— = 0 (1)
дх
Роу-
ду
дх, дхг дх,
- г ]
Росру-
дв
дх,.
■ = к
дв ду
—г + Т„—-
дх 2 - дх
(2) (3)
Здесь - компоненты вектора скорости в прямоугольном декартовом базисе, р0 - плотность (константа), р - давление, т, - компоненты тензора вязких напряжений, в - температура, ср - удельная теплоемкость при постоянном давлении, к - коэффициент теплопроводности, хг - декартовы координаты.
Рис. 1. Геометрия пресс-формы с переменным сечением для течения полимера
Дополним систему (1) определяющими соотношениями неньютоновской жидкости. Для этой цели могут быть рассмотрены различные модели, например, использованы модели А1 и Ау фойгтовских изотропных вязких сред, которые предложены в [21]. В [22] показано, что для изотропной неньютоновско-вязкой среды фойгтовского типа тензор вязких напряжений можно представить в виде квазилинейной функции от тензора скоростей деформации Я- :
Т- = = V Где ц - коэффициент вязкости
( ду ду,)
— + —-кдх- дх . ,
(4)
м(в, 12 (Я) = а(в)м0(Т 2 (П)),
(5)
зависящий от температуры в и 12 (П) - второй инварианта тензора скорости деформации [23,24]:
12 (Э ЬТ2^
(6)
Рассмотрим модель Carreau-Yasuda [1,22], в которой зависимость коэффициента вязкости ¡и0(/2(Э)) имеет степенной вид
¡0 - и
1 +
X 2 )" Г
)/ а
(7)
Здесь а - безразмерный параметр, описывающий переходную область между областью нулевой скорости сдвига и областью, где имеет место степенная зависимость коэффициента вязкости от /2 (Э) ; п - показатель, характеризующий режим разрежения сдвига, наблюдаемый при высоких скоростях сдвига для псевдопластических жидкостей (0 < п < 1), а X - время релаксации - постоянная, обратной величиной к которой является критическая относительная скорость сдвига, обозначающая начало разрежения.
Температурная зависимость коэффициента вязкости описывается функцией а(в) - коэффициентом температурного сдвиг, который подчиняется закону Аррениуса [25,26]:
а
(в)-
exp
Е
(
Я
1
1
V
в-в0 ва в0 ;
(8)
Здесь Е0 и Я - постоянная энергии активации и газовая постоянная, соответственно.
Граничные условия к системе (1)-(4) были следующими: на входе слева (рис.1) в форму подается давление, а на выходе из пресс-формы задано давление. На твердых стенках пресс-формы задаются граничные условия прилипания, т.е. вектор скорости на верхней и нижней стенках формы равен нулю. На стенках пресс-формы задается фиксированная температура. Таким образом, имеем: на стенках пресс-формы: — 0 ; в — вШ1.
на входе в пресс-форму : рпШ — р(х, у) ; в — вп
Л дв
на выходе из пресс-формы: р — 0 ; —^
дП
¡пШ '
— 0.
г
Вариационная формулировка. Для применения метода конечных элементов дадим слабую (вариационную) формулировку задачи
(1)-(4) с использованием метода Галеркина, в котором весовые функции для полей давления, температуры и скорости. (у{, р,0) аппроксимируются разложениями вида
IV1
V (х, г ) = (х>? (г ) = р^ ,
т=1 L
Р(хг) = (х)Р1(г) = УТр,
I=1
0(х, г )=± ф (х0 (г )=фТ0,
(9а) (9Ь) (9с)
,=1
Где (Рт (х) , У/т (Х) фт (х) - функции формы р = р (x), Я>м (х)}, у = {у1(х),...,ум(х)}, ф = {ф1 (х),...,фм(х)} - координатные столбцы,
составленные из функций формы, V™(г), р1 (г) 0, (г) - значения узловых функций, V = М),..., vM (г)}, р = {Р1 (г),..., Рм (г)}, 0 = {01 (г),...,0м(г)} - координатные столбцы, составленные из узловых функций.
Умножая внешним образом уравнение несжимаемости (1) на столбец функций формы у и интегрируя по области Ое, получаем вариационную формулировку уравнения несжимаемости
дх,
V = 0.
(10)
Умножая внешним образом уравнение движения (2) на столбец функций формы р, интегрируя по области Ое и применяя к тензору напряжений метод интегрирования функций по частям [21, 22], получаем вариационную формулировку уравнения движения (2)
Ро
Г рдр (О дх
дх
¥ТГ] +
др дрТ —ао
дх] дх1
V +
1о.
др дрТ /и——— аО дх, дх,
У 1
V,
Г ^ао
дх,
Р = -ГГе ■.
(11)
Г
О.
где tj = п^ (-рд^ + тт) - компоненты вектора напряжения на границе области, которые содержат давление и компоненты тензора вязких напряжений. Величины ti - полагаются заданными на границе Ге области Ое.
Вариационная формулировка уравнения энергии (3) записывается аналогичным образом
РоС
г ф(рту
дх,.
е+к
+
ГоФ
V 1
Вводя обозначение для матриц
Г дф дфТ
дх 1 дх1
{дрт у+дтТ.у
дх,. i дх 1
е = Фч^
У АО дх
(12)
о, = г , д = г , бн = г р
i JQe яу 1 •Ъе ях 1 ье ^
дрт
дх.
Гое Р~дх;
др дрт
дх , дх, 1 ,
К = -
Г рё , св(у)= Г ф(рту, Бв = Г
Зте^ , 1 ^ 1! Яг 11 Ье дх дт
ув = Гое ф
(
Яр—У +Яр—у , дх 1
дх1
V 1
Ях ,
Яр
ле дх дх1
дх
У^О , Ке = фqds,
систему вариационных уравнений (10)-(12) можно записать следующим образом:
- оТУ, = о
РоОУТУ + + Бу - 0,Р = К РоСрСеУ )е+коее = уб+кв
(13)
(14)
(15)
Изопараметрические конечные элементы. Для выбора функций формы рт(х),щт(х) фт(х) - применим изопараметрический конечный элемент, который подчиняется условию ЬББ, т. е. для аппроксимации скорости и температуры в каждом КЭ используется девятиточечная квадратичная интерполяционная функция, а для аппроксимации давления используется четырехточечная линейная интерполяционная функция (рис. 2).
ш
□
ш □
Рис. 2. Изопараметрический КЭ, удовлетворяющий условию ЬББ
Безразмерная форма стандартной интерполяционной функции четырехстороннего КЭ, удовлетворяющего условию ЬББ, выглядит следующим образом
(1 -п)п
- 2(1 )(1 -Пп
-(1 + -п)п
- 2(1 -#(1 -п2)
(1 )(1 -п2)
2(1 + #(1 -п2) -(1 -^1 + п)%п 2(1 )(1+п)п
(1+^)(1+п)п
(1 -п) (1 + ^)(1 -п) (1+Й(1+п) (1 -Й(1+п)
<=4
(1 -п)п
- 2(1 )(1 -п)п
-(1 + ^)(1 -п)п
- 2(1 -#(1 -п2) (1 )(1 -п2) 2(1 + #(1 -п2)
-(1 -^)(1 + п)%п 2(1 )(1+п)п (1+й(1+п)п
(16)
где -1 <%,п -1, здесь п - безразмерные (барицентрические) координаты. 76
Якобиева матрица преобразования между от декартовых координатам к безразмерным координатам имеет следующий вид:
дх дрк
ду=дрк
Ук
дх дфк
дп дп
<У = дук дп дп
(17)
Ук
Соотношения между интерполяционными функциями в разных системах координат имеют следующий вид:
дх = J-
дРк дФк
_ ду _ дп
(18)
где I - якобиан матрицы Якоби J =
дх ду дх ду
дп дп
В переменных п численное интегрирование осуществляется по единичному квадрату.
Итерационный алгоритм решения матричной системы уравнений Явная форма нелинейной системы алгебраических уравнений (13)-(15) имеет следующий вид
Ро
яу? яу? о
02г2т о 0 0 0
РоСр [С (VI)+С 2 (¥2 )\э +К [б2п + Бв22 2 = ¥Б + ^
V" 2БП + Б 22 Б 21 - 0 'V'
¥2 + Б12 Б„ + 2Б22 - 02 ¥2 =
Р _-- от 0 Р 0
(19)
(20)
В общем виде систему нелинейных уравнений (19) можно записать в следующем виде:
К (V )¥ = Ъ (V)
(21)
х
х
к
к
где V =
V р
, а К (V) - обобщенная матрица размера 3 х 3
К (V )=р A(V )+B(V ),
QlVT Q2VT 0~ " 2ВХХ + ^22 А1 - О
A(V) = QlV2T Q2V2T 0 , В (V) = А2 Ви + 2 Б22 -О2
0 0 0 - ОТ - ОТ 0
а ¥ (V) =
¥2 0
вектор правых частей.
Для повышения точности решения нелинейно системы (19)-(20), рассмотрим метод Ньютона-Рафсона со вторым порядком точности [19]. Для этого перепишем уравнение. (21) в виде
R(v ) = К(V V - ¥(V) = 0
(22)
Метод Ньютона основан на использовании усеченного ряда Тейлора R(v) в окрестности известного решения V" на предыдущем шаге итерации:
0 = R(у" )+-
дЯ
дV
ду+о^ )2
(23)
где ДУ = У"+1 - V".
Отбрасывая члены второго порядка и выше, получим
я(у" ) =
дЯ
дV
(v"+1 - V" )=- J (V" )(v"+1 - V")
V"
где J (V) - якибиева матрица
(24)
J (V" ) =
дЯ
дV
(25)
Решая уравнение (24) для V"+ , получаем
уп+1 = уп - J-1 (у. )Е(уп )
С учетом (21) и (22) матрица Я (У) имеет следующий вид
(26)
" Я" + 2Б„ + Б22 ^ + ^21 - ^ " у1 ' Р1 '
я(уп )= Я2 = ОУТ +Д2 + Б„ + 2^22 - ^ У2 - А
_ Яз _ - ^ - ^ 0 Р 0
(27)
Якобиева матрица (27) записывается следующим образом
J (Уп):
дЯ1 дЯ1 дЯ1
ду дУ дР
дЯ2 дЯ2 дЯ2
ду дУ2 дР
дЯ3 дЯ3 дЯ3
ду дУ2 дР
где
^ = 2&УГ + 02У2Г + 2Пи + ^22, ^ = 02^ + А1 ду дУ2
Я = ОУ! + А, ^ = е1У1г + 202У2Г + А1 + Ю2
дУ
дУ
дЯ3 дУ
дЯ
1 дУ
Я
дР
= 0.
дя_
дР
Я
дР
= -ор
= -02 Р
Для решения уравнения энергии (20) вводится параметр в корректировки скорости сходимости итерационного алгоритма
РоСР [С (VI)+Свв V )][(1 - РУ + рвп+1 ]+ фв + Бв22 ]^п+1 = УБ + ^ (28)
где Ре [0,1].
Численный алгоритм решения состоит из следующих этапов: 1) на основании заданных начальных условий рассчитывается начальная вязкость, номер итерации устанавливается равным п = 1, т.е. лР1] = V, р(п> = р0, 2(п) = 20.
2) для граничных условий к системе уравнений Навье-Стокса неньютоновской жидкости применяется итерационный алгоритм Нью-тона-Рафсона, получены у|"+1), р(п+1).
3) расчитываются скорость и давления на п-м шаге итерации по формуле (26) и подставляются в уравнение энергии (28) для нахождения распределения температуры 2(п+1).
4) осуществляется проверка сходимости на п -м шаге итерации,
если (X(п+1)-X(п)УX(п) <е, X = V.,р,2 затем итерационный про-
^ 'ад
цесс завершается; в противном случае п +1 ^ п, шаг 2) перезапускается.
Исходные данные для численного моделирования неизотермического потока полимерного раствора. Для рассматриваемой области решения задачи течения неньютоновской жидкости в пресс-форме переменного сечения (рис.1) была построения дискретная сетка из изопараметрических четырехугольных элементов. Сетка показана на рис. 3.
Рис. 3. Дискретная сетка для задачи течения неньютоновской жидкости в
пресс-форме
Значения реологических констант неньютоновской жидкости , соответствующие модели (5)-(8), приведены в табл. 1.
Температура стенки полости 2Ш1 составляет 250°С, а начальная температура материала 2М составляет 180°С.
В неизотермическом потоке характеристическое калориметрическое число Пекле Ре отражает соотношение между термическим конвекционным термом и диффузионным термом, и определяется следующим образом:
Ре = ,
К
где V и I - характерная скорость и длина соответственно, в данной задаче приняты V = 0.05 ти I = 0.2т.
В задаче проводился анализ влияния числа Пекле на результаты решения. Были проведены расчеты при 3 различных числах Пекле, которые были выбраны равными Ре = 10, 50, 250. При фиксированных параметрах задачи, три значения числа Ре соответствовали 3 значениям теплоемкости с (см. табл.1).
Таблица 1
Константы модели неньютновской жидкости [27]
Параметр Значение Параметр Значение
р0 - плотность (кг/м3) 921 /Л0 - вязкость при нулевой скорости сдвига (Рас) 10000
с - удельная теплоёмкость Дж/(кгК) 2, 10, 20, 50 - вязкость при бесконечной скорости сдвига (Ра с) 500
к - коэффициент теплопроводности Вт/(мК) 1.842 п - показатель степени в уравнении нелинейной вязкости 0.45
в0 - начальная температура °С -273 X - время релаксации 2.5
да - эталонная температура °С 205 а - безразмерный параметр 2
Я - газовая постоянная Дж^моль-К). 8.31446 Е0 - постоянная энергии активации Дж^м) 18800
Таблица 2
Параметры пресс-формы для течения полимера
Параметр Значение Параметр Значение
Высота Н (м) 0.2 Длина L (м) 0.4
Угол наклона у (°) 60 Коэффициент сужения 0.125, 0.25, 0.3, 0.5
Анализ влияние числа Пекле. Давление на входе в пресс-форму было задано р = 105Ра, коэффициент фиксированной усадки составлял 0.25. На рис. 4 показано распределение скорости, давления, вязкости и температуры жидкости при Ре = 250.
На рисунке 5 показаны соответствующие распределения параметров при Ре = 50, а на рис.6 - при числе Ре = 10.
Распределение скорости, давления, вязкости и температуры симметрично относительно центральной оси полости. Анализируя приведенные выше численные результаты, можно сделать вывод, что влияние числа Пекле на скорость потока не очень значительно: увеличение
числа Пекле усиливает влияние тепловой конвекции, а уменьшение числа Пекле усиливает влияние теплопроводности. В то же время влияние числа Пекле на вязкость очень существенно: значение скорости в самой узкой части пресс-формы является наибольшим, и изменение вязкости, зависящей от скорости, в этой области является наиболее значительным.
д
Рис. 4. Распределение:
а - горизонтальной скорости; б - вертикальной скорости; в - давления; г - вязкости; д - температуры
д
Рис. 5. Распределение:
а - горизонтальной скорости; б - вертикальной скорости; в - давления; г - вязкости; д - температуры
а
б
в
г
Рис. 6. Распределение:
а - горизонтальной скорости; б - вертикальной скорости; в - давления; г - вязкости; д - температуры
а
б
в
г
д
Анализ влияние коэффициента сужения. Было проведено исследование влияния коэффициента усадки пресс-формы на скорость потока жидкости и теплопередачу. Значение числа Pe было фиксировано на уровне 100, а внешнее давление составляло 105Pa. Коэффициент сужения изменялся и составлял в расчетах: 0.125, 0.25, 0.5.
Результаты численного для значения коэффициент сужения 0,125 представлены на рис. 7. На рис. 8 показаны соответствующие значения параметров течения жидкости при коэффициенте сужения 0,25, а на рис. 9 - для коэффициента сужения 0,5:
velocity U хЮ"
8
6
д
Рис. 7. Распределение:
горизонтальной скорости; б - вертикальной скорости; в - давления; г - вязкости; д - температуры
v3
а
Рис. 8. Распределение:
а - горизонтальной скорости; б - вертикальной скорости; в - давления; г - вязкости; д - температуры
а
б
в
г
д
velocity U
ш
я
0.05
velocity V
ЦТ ■
б
Pressure Р
001
о
-0.01
хЮ4
>—
Viscosity
Temperature
12000 10000 8000 6000 4000 2000
250
200
Рис. 9. Распределение:
а - горизонтальной скорости; б - вертикальной скорости; в - давления; г - вязкости; д - температуры
Анализ приведенных выше численных результатов показывает, что увеличение коэффициента сужения увеличивает поток жидкости. Увеличение скорости приводит к увеличению числа Пекле, поэтому распределение поля температуры показывает, что увеличение коэффициента сужения увеличивает конвекцию жидкости.
Интересным результатом является распределение коэффициента вязкости. На графике с коэффициентом сужения 0.5, в самой узкой об-
а
в
г
ласти пресс-формы происходит увеличение вязкости. Это связано с нелинейной зависимостью вязкости от температуры и скорости. Поскольку скорость течения в узкой части пресс-формы наибольшая, то и вязкость жидкости в этой зоне повышается по сравнению с широкой областью пресс-формы.
Анализ влияния внешнего давления. Наконец, был исследовано влияние внешнего давления на течение и теплопередачу в полимерном растворе. Расчеты проводились при следующих значениях параметров: коэффициент сужения 0.3, число Pe = 100, внешнее давление принимало значения: 5 х104Pa, 105Ра и 2 х105 Pa.
Результаты моделирования для случая внешнего давления 5 х 104 Ра показаны на рис. 10.
velocity U
0.01
velocity V
б
Pressure Р
х1ЕГ
Viscosity
Temperature
I
250
200
Рис. 10. Распределение:
а - горизонтальной скорости; б - вертикальной скорости; в - давления; г - вязкости; д - температуры
а
в
г
На рис. 11 и 12 показаны результаты моделирования для случаев значений внешнего давления 105 Ра и 2x105 Ра , соответственно.
velocity U
0.04 0.02
а
б
в
г
д
Рис. 11. Распределение:
горизонтальной скорости; б - вертикальной скорости; в - давления; г - вязкости; д - температуры
а
Рис. 12. Распределение:
а - горизонтальной скорости; б - вертикальной скорости; в - давления; г - вязкости; д - температуры
а
б
в
г
д
Влияние внешнего давления в основном сказывается на скорости течения жидкости: увеличение внешнего давления приводит к увеличению скорости жидкости, тем самым, вызывая увеличение конвекции жидкости.
Выводы. Разработана математическая модель неизотермического потока неньютоновских жидкостей в областях сложной формы с переменными границами. Получена вариационная формулировка задачи. Предложен численный итерационный метод решения задачи течения неньютоновских жидкостей, основанный на итерационном методе Ньютона-Рафсона. Результаты численного моделирования демонстрируют существенное влияние реологических свойств жидкости, внешнего давления и характеристик формы на неизотермический поток неньютоновских жидкостей, что необходимо учитывать при технологических процессах переработки пластмасс.
ЛИТЕРАТУРА
[1] Bird R.B, Armstrong R.C., Hassager O., Dynamics of polymeric liquids. Vol. 1: Fluid mechanics. John Wiley & Sons, 1987, 649p.
[2] Chhabra R.P., Richardson J.F., Non-Newtonian Flow: Fundamentals and Engineering Applications. Elsevier, 1999, 436p.
[3] Larson R.G., The structure and rheology of complex fluids (topics in chemical engineering). Oxford University Press, 1999, 663p.
[4] Димитриенко Ю.И., Захарова Ю.В., Богданов И.О. Математическое и численное моделирование процесса фильтрации связующего в тканевом композите при RTM методе изготовления. Университетский научный журнал, 2016. №19. С. 33-43.
[5] Димитриенко Ю.И., Шпакова Ю.В., Богданов И.О., Сборщиков С.В. Моделирование процесса многоуровневой фильтрации жидкого связующего в тканевом композите при RTM-методе изготовления. Инженерный журнал: Наука и инновации. М.: МГТУ им. Н.Э. Баумана, 2015. № 12(48). 7 с.
[6] Ю.И. Димитриенко, И.О. Богданов Многомасштабное моделирование процессов фильтрации жидкого связующего в композитных конструкциях при RTM методе изготовления. Математическое моделирование и численные методы, 2017. № 2. С. 3-27.
[7] Димитриенко Ю.И., Иванов М.Ю. Моделирование нелинейных динамических процессов переноса в пористых средах. Вестник МГТУ им.Н.Э.Бау-мана. Сер. Естественные науки, №1. 2008. С.39-56.
[8] Dimitrienko Yu.I., Dimitrienko I.D. Simulation of local transfer in periodic porous media. European Journal of Mechanics/B-Fluids, 2013. № 1. P.174-179.
[9] Димитриенко Ю.И., Левина А.И., Боженик П. Конечно-элементное моделирование локальных процессов переноса в пористых средах. Вестник МГТУ им.Н.Э.Баумана. Сер. Естественные науки, 2008. № 3.С.90-104
[10] Ю.И. Димитриенко, И.О. Богданов Многомасштабное моделирование процессов фильтрации в пористых средах. Инженерный журнал: наука и инновации, 2018. № 3(75). 19 с.
[11] Zienkiewicz O.C., Taylor R.L., Zhu J.Z., The Finite Element Method: Its Basis and Fundamentals: Its Basis and Fundamentals. Elsevier, 2005, 733p.
[12] Zienkiewicz O.C., Cheung Y.K., Finite elements in the solution of field problems. The Engineer, 1965, vol. 220, iss. 5722, pp. 507-510.
[13] Zienkiewicz O.C., Taylor R.L., Nithiarasu P., The finite elements methods for fluid Dynamics. Elsevier, 2005, 435p.
[14] Oden J.T., The finite element method in fluid mechanics. Finite element methods in continuum mechanics, 1973, pp. 151-186.
[15] Lewis R.W., Nithiarasu P,. Seetharamu K.N., Fundamentals of the finite element method for heat and fluid flow. John Wiley & Sons, 2004, 341p.
[16] Nassehi V., Practical aspects offinite element modelling ofpolymer processing. John Wiley & Sons, 2002, 273p.
[17] Han X.H., Li X.K., An iterative stabilized CNBS-CG scheme for incompressible non-isothermal non-Newtonian fluid flow. International journal of heat and mass transfer, 2007, vol. 50, iss. 5-6, pp. 847-856.
[18] Mu Y., Zhao G.Q., Wu X.H., Zhai J.Q., Modeling and simulation of three-dimensional planar contraction flow of viscoelastic fluids with PTT, Giesekus and FENE-P constitutive models. Applied Mathematics and Computation, 2012, vol. 218, iss. 17, pp. 8429-8443.
[19] Reddy J.N., Gartling D.K., The finite element method in heat transfer and fluid dynamics. CRC press, 2010, 489p.
[20] Dimitrienko Yu.I. Thermomechanics of composites under high temperatures. Springer, 2015, 434p.
[21] Димитриенко Ю.И. Механика сплошной среды. В 4 т. Т. 2. Универсальные законы механики и электродинамики сплошных сред. Москва, Изд-во МГТУ им. Н.Э. Баумана, 2011, 560 p.
[22] Димитриенко Ю.И. Нелинейная механика сплошной среды. Москва, Физма-тлит, 2009, 624p.
[23] Димитриенко Ю.И. Механика сплошной среды. В 4 т. Т. 1: Тензорный анализ. Москва, Изд-во МГТУ им. Н.Э. Баумана, 2011, 367 p.
[24] Dimitrienko Yu.I. Tensor analysis and nonlinear tensor functions. Springer, 2002, 662p.
[25] Han C.D., Rheology and processing ofpolymeric materials: Volume 1: Polymer Rheology. Oxford University Press on Demand, 2007, 707p.
[26] Peters G.W.M., Baaijens F.P.T., Modelling of non-isothermal viscoelastic flows, J. Non-Newton. FluidMech. 1997, vol. 68, pp. 205-224.
[27] Toth G., Bata A., Belina K. Determination of polymer melts flow-activation energy a function of wide range shear rate. IOP Conf. Series: Journal of Physics: Conf. Series 1045 (2018) 012040
Статья поступила 12.06.2018
Ссылку на эту статью просим оформлять следующим образом:
Димитриенко Ю.И., Шугуан Ли Конечно-элементное моделирование неизотермического стационарного течения неньютоновской жидкости в сложных областях. Математическое моделирование и численные методы, 2018, № 2, с. 70-95.
Димитриенко Юрий Иванович - д-р физ.-мат. наук, профессор, директор Научно-образовательного центра «Суперкомпьютерное инженерное моделирование и разработка программных комплексов» (НОЦ «СИМПЛЕКС») МГТУ им. Н.Э. Баумана; заведующий кафедрой «Вычислительная математика и математическая физика» МГТУ им. Н.Э. Баумана. e-mail: dimit.bmtstu@gmail.com
Шугуан Ли - аспирант кафедры «Вычислительная математика и математическая физика» МГТУ им. Н.Э. Баумана. e-mail: shuguangli2008@gmail.com
Mathematical simulation of non-isothermal steady flow of non-Newtonian fluid by finite element method
© Y.I. Dimitrienko, Shuguang Li
Bauman Moscow State Technical University, Moscow, 105005, Russia
The finite element method is used to simulate the nonisothermal flow of non-Newtonian viscous fluids in complex geometries. The Carreau-Yasuda model of a non-Newtonian fluid is considered, in which the dependence of the viscosity coefficient on the second invariant of the strain rate tensor has a power form. A variational formulation of the problem of the motion of a non-Newtonian fluid for a plane case is obtained. The iteration algorithm of Newton-Raphson is used to solve the Navier-Stokes equations system, and the Picard iteration algorithm is used to solve the energy equation. The problem of the movement of a polymer mass in a mold of complex variable cross section in the presence of an uneven temperature field is considered. With the help of finite element modeling, a numerical analysis of the effect of various parameters on the movement of a liquid and the heat transfer of a polymer material at different values of external pressure was carried out. It is shown that the nature of the motion of a non-Newtonian fluid essentially depends on the rheolog-ical properties of the fluid and the characteristics of the geometric shape, which must be taken into account in technological processes of plastics processing.
Keywords: finite element method; non-newtonian fluid; non-isothermal flow; iterative algorithm; complex geometric shapes, convection.
REFERENCES
[1] Bird R.B, Armstrong R.C., Hassager O., Dynamics of polymeric liquids. Vol. 1: Fluid mechanics. John Wiley & Sons, 1987, 649p.
[2] Chhabra R.P., Richardson J.F., Non-Newtonian Flow: Fundamentals and Engineering Applications. Elsevier, 1999, 436p.
[3] Larson R.G., The structure and rheology of complex fluids (topics in chemical engineering). Oxford University Press, 1999, 663p.
[4] Dimitrienko Yu.I., Zakharova Yu.V., Bogdanov I.O. Mathematical and numerical simulation of the binder filtration process in the fabric composite under the RTM manufacturing method Humanities and Science University Journal, 2016. no. 19. PP. 33-43.
[5] Dimitrienko Yu.I., Shpakova Yu.V., Bogdanov I.O., Sborshchikov S.V. Inzhe-nernyj zhurnal: Nauka i innovacii - Engineering Journal: Science and Innovation, 2015, no. 12(48), 7 p.
[6] Dimitrienko, I.O. Bogdanov Multiscale modeling of filtration liquid binding processes in composite designs at RTM production method. Mathematical Modeling and Computational Methods, 2017. no. 2. pp. 3-27.
[7] Dimitrienko Yu.I., Ivanov M.Yu. Modeling of Nonlinear Dynamical Processes of Transfer in Porous Media. VestnikMGTU im N.E.Baumana. Ser. Estestvennye nauki - Herald of the Bauman Moscow State Technical University. Natural Sciences,, 2008. no. 1. PP.39-56.
[8] Dimitrienko Yu.I., Dimitrienko I.D. Simulation of local transfer in periodic porous media European Journal of Mechanics/B-Fluids, 2013. № 1. P.174-179.
[9] Dimitrienko Yu.I., Levina A.I., Bozhenik P. Vestnik MGTU im N.E.Baumana. Ser. Estestvennye nauki - Herald of the Bauman Moscow State Technical University. Natural Sciences, 2008. no. 3. PP.90-104
[10] Dimitrienko Yu. I., Bogdanov I.O. Inzhenernyy zhurnal: nauka i innovatsii - Engineering Journal: Science and Innovation, 2018, no. 3(75).
[11] Zienkiewicz O.C., Taylor R.L., Zhu J.Z., The Finite Element Method: Its Basis and Fundamentals: Its Basis and Fundamentals. Elsevier, 2005, 733p.
[12] Zienkiewicz O.C., Cheung Y.K., Finite elements in the solution of field problems. The Engineer, 1965, vol. 220, iss. 5722, pp. 507-510.
[13] Zienkiewicz O.C., Taylor R.L., Nithiarasu P., The finite elements methods for fluid Dynamics. Elsevier, 2005, 435p.
[14] Oden J.T., The finite element method in fluid mechanics. Finite element methods in continuum mechanics, 1973, pp. 151-186.
[15] Lewis R.W., Nithiarasu P,. Seetharamu K.N., Fundamentals of the finite element method for heat and fluid flow. John Wiley & Sons, 2004, 341p.
[16] Nassehi V., Practical aspects of finite element modelling of polymer processing. John Wiley & Sons, 2002, 273p.
[17] Han X.H., Li X.K., An iterative stabilized CNBS-CG scheme for incompressible non-isothermal non-Newtonian fluid flow. International journal of heat and mass transfer, 2007, vol. 50, iss. 5-6, pp. 847-856.
[18] Mu Y., Zhao G.Q., Wu X.H., Zhai J.Q., Modeling and simulation of three-dimensional planar contraction flow of viscoelastic fluids with PTT, Giesekus and FENE-P constitutive models. Applied Mathematics and Computation, 2012, vol. 218, iss. 17, pp. 8429-8443.
[19] Reddy J.N., Gartling D.K., The finite element method in heat transfer and fluid dynamics. CRC press, 2010, 489p.
[20] Dimitrienko Yu.I. Thermomechanics of composites under high temperatures. Springer, 2015, 434p.
[21] Dimitrienko Yu.I. Mekhanika sploshnoy sredy. Tom 2. Universalnye zakony mekhaniki i elektrodinamiki sploshnoy sredy [Continuum Mechanics. Vol. 2. Universal laws of mechanics and electrodynamics of continuous media]. Moscow, BMSTU Publ., 2011, 560 p.
[22] Dimitrienko Yu.I. Nelineinaya mekhanika sploshnoi sredy [Nonlinear Continuum Mechanics]. Moscow, Fizmatlit Publ., 2009, 624 p.
[23] Dimitrienko Yu.I. Mekhanika sploshnoy sredy. T.1. Tenzornoe ischislenie. [Continuum Mechanics. Vol. 1. Calculus of Tensors.]. Moscow, BMSTU Publ., 2011, 463 p.
[24] Dimitrienko Yu.I. Tensor analysis and nonlinear tensor functions. Springer, 2002, 662p.
[25] Han C.D., Rheology and processing ofpolymeric materials: Volume 1: Polymer Rheology. Oxford University Press on Demand, 2007, 707p.
[26] Peters G.W.M., Baaijens F.P.T., Modelling of non-isothermal viscoelastic flows, J. Non-Newton. FluidMech. 1997, vol. 68, pp. 205-224.
[27] Toth G., Bata A., Belina K. Determination of polymer melts flow-activation energy a function of wide range shear rate. IOP Conf. Series: Journal of Physics: Conf. Series 1045 (2018) 012040
Dimitrienko Y.I., Dr. Sc. (Phys.- Math.), Professor, Director of Research and Education Center of Supercomputer engineering modeling and program software (Simplex), Bauman
Moscow State Technical University, Head of the Department of Computational Mathematics and Mathematical Physics, BMSTU. e-mail: dimit.bmtstu@gmail.com
Shuguang Li, post-graduate student of the Department of Computational Mathematics and Mathematical Physics, BMSTU. e-mail: shuguangli2008@gmail.com
