УДК 579.64
Структурно-феноменологическая реологическая модель для инженерных расчетов течений полимерных сред
А.А. Лаас1, Г.О. Рудаков2, Г.В. Пышнограй1,3, К.Б. Кошелев4
'Алтайский государственный университет (Барнаул, Россия) 2Алтайский государственный педагогический университет (Барнаул, Россия) 3Алтайский государственный технический университет им. И.И. Ползунова (Барнаул, Россия)
4Институт водных и экологических проблем СО РАН (Барнаул, Россия)
Structural-Phenomenological Rheological Model for Engineering Calculations of Polymeric Media Flows
A.A. Laas1, G.O. Rudakov2, G.V. Pyshnograiu, K.B. Koshelev4
'Altai State University (Barnaul, Russia)
2Altai State Pedagogical University (Barnaul, Russia)
3Polzunov Altai State Technical University (Barnaul, Russia)
4Institute for Water and Environmental Problems SB RAS (Barnaul, Russia)
В области нелинейной вязкоупругости изучение поведения течений растворов и расплавов полимеров позволяет более подробно описать реологические свойства и точнее оценить адекватность реологических моделей. Для описания реологического поведения расплавов разветвленных полимеров предложена новая структурно-феноменологическая модель, которую можно рекомендовать для инженерных расчетов течений полимерных сред. На основе модели, которая получена, исходя из модифицированной модели Виноградова-Покровского, основанной на микроструктурном подходе к описанию динамики полимерной жидкости, были рассчитаны стационарные вискозиметрические функции при простом сдвиге и одноосном растяжении: стационарная сдвиговая вязкость, коэффициент первой разности нормальных напряжений, стационарная вязкость при одноосном растяжении. Также было исследовано влияние параметров модели на вид этих зависимостей. Показано, что модель с хорошей точностью описывает нелинейное вязкоупругое поведение текучих полимерных систем: аномалию вязкости, падение коэффициента первой разности нормальных напряжений, немонотонный характер зависимости стационарной вязкости при растяжении от скорости растяжения. Проведено сравнение вискозиметрических функций с экспериментальными данными для расплава промышленного образца полиэтилена.
Ключевые слова: реология, полимеры, расплавы полимеров, полимерная жидкость, реологическая модель, мезоскопический подход, вискозиметрические функции, сдвиговая вязкость, вязкоупругость, нелинейные эффекты, одноосное растяжение, вязкость при растяжении, щелевой канал.
DOI: 10.14258/izvasu(2022)4-18
Studying the behavior of flows of solutions and polymer melts in the field of nonlinear viscoelasticity allows describing the rheological properties in more details and more accurately assess the adequacy of rheological models. A new structural-phenomenological model is proposed to describe the rheological behavior of melts of branched polymers. This model can be recommended for engineering calculations of flows of polymeric media. The model is obtained from the modified Vinogradov-Pokrovsky model which is based on the microstructural approach and describes the dynamics of a polymer fluid. Stationary viscometric functions for simple shear and uniaxial tension, as well as stationary shear viscosity, the first-difference coefficient of normal stresses, stationary viscosity at uniaxial tension, are calculated using the obtained model. Also, the influence of model parameters on the form of the functions has been studied. It is shown that the model describes with good accuracy the nonlinear viscoelastic behavior of flowing polymer systems: an anomaly in viscosity, a drop in the coefficient of the first difference of normal stresses, and a nonmonotonic nature of the dependence of the steady-state elongation viscosity on the tensile rate. The viscometric functions data are compared with the experimental data for an industrial polyethylene melt sample.
Key words: rheology, polymers, polymer melts, polymer fluid, rheological model, mesoscopic approach, viscometric functions, shear viscosity, viscoelasticity, nonlinear effects,
uniaxial elongation, tensile viscosity, slot channel.
Введение
В настоящее время в области моделирования технологических процессов переработки структурно-неоднородных жидкостей, к которым относятся растворы и расплавы полимеров, нефть, краски и др., сложилась противоречивая ситуация. Хотя эти среды и обладают рядом уникальных свойств, таких как вязкоупругость, различные формы нелинейного поведения, для описания которых уже сформулированы пригодные реологические модели в случаях простых с гидродинамической точки зрения, в инженерных расчетах в основном используют устаревшие модели [1]. К таким моделям относятся: Освальда де Валле, Карро, Карро-Ясуда, Кросса, Кросса-Вильямсона, Эллиса, модель вязкопластич-ной жидкости Бингама и др.
Применение математического моделирования при проектировании и производстве изделий из полимерных материалов определяется реологической моделью, которая должна быть достаточно простой, но отражать все соответствующие характеристики процесса [2-7].
В настоящее время предложено большое число реологических моделей, учитывающих сложное строение полимерных молекул, однако в литературе упоминается несколько реологических моделей, которые с точностью описывают нелинейные эффекты при одноосном растяжении и простом сдвиге, поэтому возникает потребность в проведении и последующей интерпретации реологических измерений с более сложной структурой течений. Вместе с тем к инженерным моделям предъявляют требования простоты и надежности в вычислениях [8, 9].
Все модели, описывающие реологическое поведение полимерных жидкостей, разделяют на два принципиально разных класса, которые отличаются подходами к описанию процессов, — это феноменологический и микроструктурный подходы. В фено-
менологическом подходе динамика макроскопических тел строится на основе общих закономерностей, которые обнаруживаются из опыта. К феноменологическому типу относят модели Максвелла [10], Олдройда [11], Прокунина-Леонова [5, 12], а также модель К-ВК2 [13].
Другой класс моделей строится на основе мезоско-пического подхода [13-15]. В таких моделях динамика макромолекул описывается на основе модельных представлений и поэтому учитывается некоторое приближение, например, строение молекул полимера и процессы межмолекулярного взаимодействия. Из реологических моделей мезоскопического подхода чаще всего используются модель РТТ [7], модель Гизекуса [4], модель Рот-Рот [6] и Виноградова-Покровского [2, 13-15].
Большинство существующих реологических моделей дает хорошие результаты для описания виско-зиметрических функций, поэтому их можно использовать для моделирования более сложных течений растворов и расплавов полимеров. Однако проведение расчетов в областях с нерегулярной геометрией может приводить к численной неустойчивости, что требует применения достаточно сложных методов регуляризации. При использовании таких моделей в инженерной практике можно столкнуться с ситуацией, когда незначительное изменение параметров течения приведет к невозможности получить решение.
Реологическая модель для инженерных расчетов
Сформулируем реологическую модель, пригодную для использования в инженерных расчетах, опираясь на тУР модель, которая имеет вид [2, 13-15]:
а = - р1 + ^
т
а-1№ а■ I 3
(1)
й „ _ ,т 1 + (к - в)£т а —а - Уу ■ а - а(Уу)7 +--—-—-
й тп
2 ,в а = — 7-3—а ■ а . 3т
(2)
Здесь а является тензором полных напряжений; р — гидростатическое давление; I — единичный тензор; п0 и т0 — начальные значения сдвиговой вязкости и времени релаксации; у — скорость; у — симметричный тензор скоростей деформации; оператор 4 означает взятие градиента; а — симметричный тензор наведенной анизотропии; & а — след безразмерного тензора а; к и в — параметры наведенной анизотропии, которые учитывают в уравнениях динамики макромолекулы размеры и форму макромолекуляр-ного клубка соответственно.
Одним из недостатков этой модели является наличие множественности стационарных решений,
что связано с наличием в уравнении (2) слагаемых а • а и У у ■ а. При проведении расчетов течений в областях со сложной геометрией это может приводить к неустойчивости численных алгоритмов. Попробуем сформулировать, опираясь на систему уравнений (1), (2), новую реологическую модель, которая будет лишена такого недостатка. Тензор а зависит от Уу и t. Выделив стационарную часть, представим а в виде:
а(Уу^) = в(Уу) + е(Уу^ ), где
(3)
1 + (к-в)к 5 2 ^ _ ,Г „в
—^-^-5 = - 7 + ^ V • 5 + у)7 -3 — 5 • 5
т„ 3 т„
(4)
й ^ чГ 1 + (к- 0)*г 5 „в/
—е-^У • 5 - 5 • (Уу)г +——-—-е = -3—(а • е + е • а).
й тп тп
(5)
В (5) опущены слагаемые, пропорциональные е • е. Перепишем (4) в виде
1+к о)
5 = -7 + т0 • 0ДО)(ш • 7 + 7• ш ) + т0 • 02(0)7• 7, т3
(6)
где О = (т0)2 & (у • у) — это безразмерный первый ин-
Таким образом, система уравнений (1), (3), (5), (6)
вариант квадрата тензора у. Выражение (6) является представляет собой реологическую модель, которая
удобной аппроксимацией зависимостей, порождае- удовлетворяет принципу материальной объективно-
мых уравнением (4), которое может иметь несколь- сти Олдройда [3]. Если продолжить процедуру упро-
ко решений [16]. щения, то вместо (5) и (6) можно записать:
й 1 + (к О)п 2
а +--0-а = -7 + т0 • о1(О)(ш • 7 + 7• шГ) + т0 • о2(О)7• 7,
Л
(7)
где к0, п0 , сДО), с2(О) — параметры модели. к0, п0 — константы, с1(О), с2(О) — функции инварианта О, которые определим далее.
Система уравнений (1), (7) уже не является реологическим определяющим соотношением в классическом понимании, так как для нее принцип материальной объективности [3] не выполняется. Однако ее стационарные решения совпадают со стационарными решениями системы (1), (3), (5), (6) и погрешность, допущенная при замене (5), (6) на (7), будет проявляться только на начальном участке временной зависимости компонент тензора а и будет пропорциональна О или квадрату градиента скорости. Поэтому уравнение (7) будет представлять собой аналог рассматриваемой системы. Следует ожидать, что модель (7) будет успешно использована при проведении инженерных расчетов течений полимерных жидкостей.
Аналитическое и численное исследование
реологической модели
Рассмотрим реологический смысл параметров модели с1(О), с2(О). Определим их так, чтобы рассчитанные на ее основе стационарные (в этом случае а(уУ) = э(Уу) вискозиметрические функции адекватно описывали поведение реальных полимерных сред.
Для этого сначала рассмотрим простое сдвиговое течение. В этом случае только одна компонента тензора градиентов скорости V = е отлична от нуля. Опуская детали расчетов, которые приведены в работе [16], получим:
г,(О) =
П0
1 + (К0 • О)п
01(О)=2• ^^ АО» = 2П0 • т0
ф2(О) = 3П0 • т,
3 1 + (к1 • О)п Г1~ ' 1 + (к1 • О)п 01 • (О) 1 02 • (О)
1 + (К0 • О)п0 2 1 + (К0 • О)п
Если рассмотреть на основе модели (1), (7) одноосное растяжение, то, следуя работе [16], получим:
А(е) = 3^0
1
3 02 • (О) • (т0 • 5)
1 + (К0 • О)п 4 1 + (К0 • О)п
Введем безразмерную скорость растяжения 5 = т0 •е;
А( 5) = 3 П0
1
3 02 • (О) • 5
1 + (К0 • О)п 4 1 + (К0 • О)п
Ч
1 + (К0 •-52Т"
(1 + (1-в)
1 + (К0 •-52Г2
Таким образом, на основе (7) модели можно рас- следовать влияние параметров модели на вид гради-считать стационарные вискозиметрические функции ентных зависимостей: при простом сдвиге и одноосном растяжении и ис- — стационарной сдвиговой вязкости:
0
5
ф) =
По
1 + (Ко ■ 52/2Г
(8)
— коэффициента первой разности нормальных напряжений:
№) =
2Пот о
1 + (к1 ■ 52/2)я
(9)
— стационарной вязкости при одноосном растяжении:
А« = -
3По
1 + (Ко ■-5Т
1 + -
(1-в) ■ 5
1 + (К2 ^52)"
(10)
где 5 = т0у — безразмерная скорость сдвига или
5 = т0£ — безразмерная скорость растяжения. Из выражений (8)-(10) видно, что реологическое поведение системы характеризуется параметрами к0, п , к , п1, к2 , п2 и в. Также из (8)-(10) можно сделать вывод, что п и у — убывающие функции скорости сдвига, X демонстрирует немонотонную зависимость; сначала вязкость при растяжении является возрастающей функцией скорости растяжения, а затем после перехода через максимум убывает.
Рассмотрим, как влияют параметры к0, п0, к1 , п1, к2 , п2 и в на вид вискозиметрических функций (9)-(11). Как видно из выражения (8), параметры к0, п0 влияют только на стационарную сдвиговую вязкость. Характер их влияния виден на рисунках 1 и 2, откуда следует, что при малых 5 п(5) принимает постоянные значения — п0, с ростом 5 п($) убывает с постоянным наклоном. При увеличении к0 зависимости п($) сдвигаются влево без изменения своего наклона. С ростом п0 (рис. 2) наклон кривых увеличивается и при этом сдвига зависимостей не наблюдается.
ю'3 1|},! 1а"' ш* 'о шг ю' ю4
v.
Рис. 1. Влияние параметра модели к0 на зависимость сдвиговой вязкости от градиента скорости
Рис. 2. Влияние параметра модели п на зависимость сдвиговой вязкости от градиента скорости
Как видно из выражения (9), параметры к1, п1 влияют только на стационарный коэффициент первой разности нормальных напряжений у1. Характер их влияния виден на рисунках 3 и 4, откуда следует, что он подобен влиянию параметров к0, п0 на сдвиговую вязкость (рис. 1, 2). То есть при малых 5 у1(5)
принимает постоянные значения — 2 п0т0, с ростом 5 У1(5) убывает с постоянным наклоном. При увеличении к1 зависимости у1(5) сдвигаются влево без изменения своего наклона. Как видно из рисунка 4, параметр п1 ведет себя аналогично параметру п0.
Рис. 3. Влияние параметра модели к1 на зависимость коэффициента первой разности нормальных напряжений от градиента скорости
Рис. 4. Влияние параметра модели п1 на зависимость коэффициента первой разности нормальных напряжений от градиента скорости
На рисунках 5, 6 и 7 показано влияние параметров к2, п2 и в на вид зависимости стационарной элонгаци-онной вязкости Х(5) от скорости растяжения 5. Из рисунка 5 видно, что при уменьшении к2 на зависимости Х(5) появляется характерный максимум, который наблюдается в экспериментах для расплавов линейных полимеров. С ростом к2 этот максимум сглаживается и может совсем исчезнуть. Такая картина характер-
на для расплавов разветвленных полимеров. Влияние параметра п2, которое можно наблюдать на рисунке 6, тоже значительно. При малых п2 наблюдается рост Х(5), который с увеличением п2 переходит в убывание. На рисунке 7 демонстрируется влияние параметра в на зависимость Х(5). Видно, что этот параметр позволяет управлять величиной максимума на зависимости Х(5).
Рис. 5. Влияние параметра модели к2 на зависимость вязкости при растяжении от безразмерной скорости растяжения
Рис. 6. Влияние параметра модели п на зависимость вязкости при растяжении от безразмерной скорости растяжения
Рис. 7. Влияние параметра модели в на зависимость вязкости при растяжении от безразмерной скорости растяжения
Из приведенных рисунков видно, что ц(5) — убывающая функция скорости сдвига, f1(5) также является убывающей функцией скорости сдвига, Л(5) демонстрирует немонотонную зависимость; сначала вязкость при растяжении является возрастающей функцией скорости растяжения, а затем после перехода через максимум убывает. Также из этих рисунков можно сделать вывод о том, что подбор параметров модели к0, п2,к1, п1, к2, п2 и в легко осуществить путем сопоставления расчетных зависимостей сдвиговой вязкости, коэффициента первой разности нормальных напряжений и элонгационной вязкости с экспериментальными данными. Причем параметры к0, п0 определяются при сравнении сдвиговой вязкости. Параметры к1 и п1 при сравнении с коэффициентом первой разности нормальных напряжений, а параметры к2 п2 и в — при сравнении с элонгационной вязкостью.
Проведем теперь сравнение с экспериментальными данными для промышленного образца полиэтилена Вга1еп [9], результаты которого приведены на рисунке 8. На этом рисунке точками обозначены экспериментальные данные, черными кривыми обозначены результаты расчета по новой модели (8-10), красные кривые относятся к расчетам по шУР модели (2).
Продемонстрируем преимущества использования реологической модели (1), (7) при расчетах пу-азейлевских течений в плоских каналах. На основе модели (1), (7) были выполнены расчеты напорного течения полимерной жидкости в каналах с параллельными стенками под действием постоянного
„ йр п
перепада давления А = —. В этом случае систему йг
уравнений необходимо дополнить уравнением сохранения импульса и уравнением неразрывности. Три
Рис. 8. Сравнение экспериментальных и теоретических зависимостей сдвиговой вязкости, вязкости при растяжении и коэффициента первой разности нормальных напряжений
компоненты вектора скорости поле давления и шесть компонент тензора дополнительных напряжений при этом являются функциями двух пространственных координат. Для сравнения были взяты канал с квадратным сечением (14x14) и щелевой канал с сечением (1x14). В расчетах были использованы следующие значения параметров: к0 = 0,1; п0 = 0,35; к1 = 1; п1 = 0,7; к2 = 0,2 п2 = 0,5; п0 = 5000: т0 = 0,1 [16]. Эти значения соответствуют течению расплава полиэтилена при температуре 200 "С [17].
Расчеты производились методом конечных элементов. Нелинейный характер задачи преодолевался методом последовательных исключений и разбиением на физические процессы. Поля скоростей и давлений определяли в рамках гидродинамической задачи, а в рамках упругой задачи определялись поля дополнительных напряжений. Затем итерации повторялись. При этом на каждом временном шаге для до-
стижения необходимой точности требовалось от двух до шести итераций.
Для оценки аппроксимационных свойств расчетной схемы проводились расчеты на сетках с различным числом элементов, для случая канала (14x14) брали: 590, 2360, 9400, 37760, а для щелевого канала (1x14) брали: 240, 960, 3840 и 15360. В качестве контролируемых параметров выбирали значение объемного
расхода Q = ^^ ш(х,у)йхйу и среднеквадра-
тическую скорость вихревого течения
Я = ^^ -^ы(х,у)2 + у(х, у)2 йхйу. Результаты расчетов
при А=7340 для случая (14x14) и А=385000 приведены в таблице. Эти значения были подобраны, чтобы обеспечить одинаковый расход через сечение канала.
Влияние числа расчетных ячеек на значения Q и Я
(14x14) (1x14)
число ячеек а Я число ячеек а Я
590 6683.07 1.448 240 6595 1.36
2360 6683.23 1.449 960 6689 1.46
9400 6683.23 1.449 3840 6689 1.47
37760 6683.45 1.579 15360 6689 1.47
Как видно из таблицы, расчеты на грубой сетке несущественно отличаются от случаев с большим числом элементов. Здесь существенным фактором выступает время расчета. Как при переходе от сетки, содержащей 2360 элементов, к сетке с 9400 элементами, так и от сетки с 960 к 3840 элементами расчетное время увеличивалось на два порядка. Поэтому для получения результатов вначале проводились расчеты на сетке с 2360 элементами, которые затем уточнялись на сетке с 9400 элементами для случая квадратного сечения канала и на сетке с 960 элементами, которые затем уточнялись на сетке с 3840 элементами для щелевого канала.
Устойчивость расчетов гарантируется как применением неявной схемы, так и использованием релаксации при переходе к следующему шагу итерации.
Сходимость подтверждается проведением расчетов нестационарной задачи для безразмерных времен от 0 до 100 и относительной точностью 0,001 на каждом шаге по времени.
Заключительные замечания
Полученная модель позволяет достаточно точно описывать стационарные и нестационарные характеристики расплавов разветвленных и линейных полимеров. При этом следует ожидать, что полученная здесь модель окажется пригодной и для концентрированных растворов полимеров. В дальнейшем предполагается использовать эту модель для проведения более сложных расчетов, а также при разработке численных методов трехмерных течений в качестве простой прогнозной модели для стабилизации вычислительных процедур.
Библиографический список
1. Борзенко Е.И., Хегай Е.И. Численное моделирование стационарного течения жидкости Балкли-Гершеля в канале с внезапным расширением // Вестник Том. гос. унта. Серия: Математика и механика. 2016. № 39. https://doi. org/10.17223/19988621/39/8.
2. Pokrovskii V.N. The Mesoscopic Theory of Polymer Dynamics // 2nd Edition, Springer, Berlin, 2010. https://doi. org/10.1007/978-90-481-2231-8.
3. Малкин А.Я. Исаев А.И. Реология: концепция, методы, приложения / пер. с англ. (Rheology: Concepts, Methods, and Applications). СПб., 2007. https://doi.org/10.1016/C2011-0-04626-4.
4. Giesekus H.A. Simple constitutive equation for polymer fluids based on the concept of deformation-dependent tensorial mobility // J. Non-Newtonian Fluid Mech. 1982. № 11. https:// doi.org/10.1016/0377-0257(82)85016-7.
5. Leonov A.I., Prokunin A.N. Analysis of simple constitutive equations for viscoelastic liquids // J. Non-Newton. Fluid Mech. 1992. № 42. https://doi.org/10.1007/978-94-011-1258-1_3.
6. Inkson N.J., McLeish T.C.B., Harlen O.G., Groves D.J. Predicting low density polyethylene melt rheology in elon-gational and shear flows with "pom-pom" constitutive equations // J. Rheol. 1999. № 43. https://doi.org/10.1122/L551036.
7. Phan-Thien N., Tanner R.I. A new constitutive equation derived from network theory // J. Non-Newtonian Fluid Mech. 1977. № 2. https://doi.org/10.1016/0377-0257(77)80021-9.
8. Малкин А.Я., Куличихин В.Г. Применение метода высокоамплитудных гармонических воздействий для анализа свойств полимерных материалов в нелинейной области механического поведения // Высокомолекулярные соединения. Серия А. 2014. № 56 (1). https://doi.org/10.7868/ S2308112014010039.
9. Pivokonsky R., Filip P., Zelenkova J. Two Ways to Examine Differential Constitutive Equations: Initiated on Steady or Initiated on Unsteady (LAOS) Shear Characteristics // Polymers. 2017. № 9. https://doi.org/10.3390/polym9060205.
10. Maxwell J. C. On the dynamical theory of gases // Trans. Roy. Soc. 1867. № 157. https://doi. org/10.1142/9781848161337_0014.
11. Oldroyd J.G. On the Formulation of Rheological Equations of State // Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences. 1950. № 200. https://doi.org/10.1098/rspa.1950.0035.
12. Leonov A.I., Prokunin A.N. Nonlinear Phenomena in Flows of Viscoelastic Polymer Fluids // Chapman and Hall. New York, 1994. https://doi.org/10.1007/978-94-011-1258-1.
13. Pyshnograi G.V., Gusev A.S., Pokrovskii V.N. Constitutive equations for weakly entangled linear polymers // Journal of Non-Newtonian Fluid Mechanics. 2009. № 164 (1-3). https://doi.org/10.1016/j.jnnfm.2009.07.003.
14. Мерзликина Д.А., Пышнограй Г.В., Пивоконский Р., Филип П. Реологическая модель для описания вискозиме-трических течений расплавов разветвленных полимеров // Инженерно-физический журнал. 2016. № 89 (3). https://doi. org/10.1007/s10891-016-1423-7.
15. Гусев А.С., Макарова М.А., Пышнограй Г.В. Мезо-скопическое уравнение состояния полимерных сред и описание динамических характеристик на его основе // Инженерно-физический журнал. 2005. № 78 (5). https://doi. org/10.1007/s10891-006-0009-1.
16. Лаас А.А., Пышнограй Г.В., Рудаков Г.О. Исследование реологических свойств наложения и трехмерного течения полимерного расплава в сходящемся канале на основе обобщения модели Максвелла // Механика композиционных материалов и конструкций, сложных и гетерогенных сред. 2021. № 11.
17. Hertel D., Münstedt H. Dependence of the secondary flow of a low-density polyethylene on processing parameters as investigated by laser-Doppler velocimetry // Journal of Non-Newtonian Fluid Mechanics. 2008. № 153. https://doi. org/10.1016/j.jnnfm.2007.12.004.