CC BY
8УДК 519.633.2
ЗАДАЧА ОБ ОСТАНОВКЕ ТЕЧЕНИЯ ВЯЗКОПЛАСТИЧЕСКОЙ СРЕДЫ
В КАНАЛЕ
Е. А. Муравлева
Задаче о течении вязкопластической среды в канале, ставшей уже классической, посвящено множество работ (см. обзор в [1, 2]). Отметим наиболее известные. В работе [3] предложена вариационная постановка задачи о течении в трубах произвольного поперечного сечения, доказаны теоремы существования и единственности решения и качественно исследованы его свойства. Применение вариационных неравенств к данной задаче рассматривается в [4]. Работа [5] посвящена исследованию течений более сложной вязкоупругопластической среды. Все это послужило основой для последующих публикаций.
Важной качественной особенностью задач о нестационарном движении вязкопластических сред является конечность промежутка времени затухания движения в отсутствие внешних сил, что принципиально отличается от соответствующего течения вязкой жидкости, которое затухает экспоненциально за бесконечно большое время. Рассмотрим изотермическое нестационарное ламинарное течение несжимаемой вязкопластической среды в бесконечно длинной цилиндрической трубе с поперечным сечением О и границей Г под действием перепада давления С. Необходимо найти скорость V, удовлетворяющую системе уравнений
дv 2 — ( Vv
"et-"v,'-T°v4^ij=c » i!x(0'T); (к
v = 0 на Г х (0,T).
Отметим, что первое уравнение системы (1) не имеет смысла в жесткой зоне Qo = {x\x Е Q, Vv(x) = 0}.
Таким образом, в задачах о течении вязкопластической среды характерной особенностью является необходимость строить решения в областях с неизвестной границей. Это обстоятельство создает большие трудности при построении эффективных численных методов их исследования.
В работе [3] предложено принять в качестве точной формулировки задачи (1) соответствующую вариационную постановку:
' v(t) Е H на (0,T), v(0) = vo Уп Е Щ-,
р J — v)dx + ц J Vv • V(-u — v) dx + то (^ j | V«| dx — J | Vt>| dx^j ^ С J(u — v) dx. ^
П П П П П
Неравенство в системе (2) автоматически включает в себя задачу о предельной поверхности. Эта поверхность отделяет область, где течение описывается первым уравнением (1), от области, где среда движется как твердое тело. На основе вариационной постановки (2) была доказана
Теорема [1]. Пусть Т = +оо и С = const, С < , где
f \Vu\ dx
7= inf -. (3)
«ея01(п)\{о} |N|L2(n)
Тогда если vo Е L2 (Q), то v(t) =0 для всех t, таких, что
t > Тс = -4" Ь Mo
Мо|Ы|£2(П)
7T0-CIQIV2J
(4)
где Ао — минимальное собственное значение оператора —А о однородными граничными условиями Дирихле (Ао > 0). Кроме того, выполняется асимптотическая оценка
\\ъ^)\\Ь2(п) (||г;о||ь2(п) + ехр(—у^)' (5)
где к = (Ао^)-1 (го7 — С\П\1/2).
Оценки (4), (5) не очень информативны, так как значение константы 7 для конкретных областей неизвестно. Речь идет только о качественном анализе поведения среды. Численное решение задачи о нестационарных течениях вязкопластической среды представляет большой интерес. Ввиду недифференцируемости определяющих соотношений она является достаточно сложной. Необходимо провести предварительный анализ предлагаемых разностных схем для решения соответствующей вспомогательной стационарной подзадачи. Следует отметить, что в недавних работах численное моделирование проводилось лишь для труб с круглым и кольцевым поперечным сечением. Это более простой случай, так как в связи с осевой симметрией задача становится одномерной. Кроме того, отсутствуют некоторые характерные эффекты течения бингамовской жидкости (в частности, появление застойных зон).
Основная сложность при численном моделировании течения вязкопластической среды обусловлена сингулярностью определяющих соотношений и невозможностью определить напряжения в тех областях, где скорость деформации равна нулю. Существуют две основные группы методов для преодоления математических сложностей, связанных с течением вязкопластической жидкости. Первая — методы регуляризации, широко используемые на протяжении многих лет. Регуляризация заключается в аппроксимации недифференцируемых определяющих соотношений гладкой функцией и в дальнейшем решении задачи о течении жидкости с нелинейной "эффективной" вязкостью. Для регуляризованных моделей не определено понятие жесткой зоны, и она вводится условием малости скоростей деформаций или критерием Мизеса для напряжений: \т\ = то. Это иногда приводит к неточному решению. Недостатки регуляризованных моделей становятся наиболее очевидными для задач об остановке течения, так как они дают качественно неверный результат: течение продолжается бесконечное время с маленькой, но ненулевой скоростью. Разумеется, при использовании указанных моделей остановкой течения можно считать достижение скоростью некоторого предельного значения. Тем не менее для данной задачи регуляризованные модели воспроизводят скорее поведение вязкой жидкости, а не вязкопластической среды.
Альтернативой могут служить вариационные методы, впервые предложенные для вязкопластиче-ских сред Ильюшиным [6]. С их помощью в работе [3] проведено исследование качественных особенностей течений в цилиндрических трубах. Применение теории вариационных неравенств к задачам течения жидкостей Бингама изложено в книге [4]. В монографии [7] предложены численные методы решения вариационных неравенств для среды Бингама (нерегуляризованная модель), основанные на применении множителей Лагранжа. Однако они казались слишком сложными и не получили должного распространения. В последние годы к вариационным неравенствам прибегают все чаще [1]. Отметим, что при использовании вариационной модели бингамовских сред для численных расчетов методом дискретизации традиционно служит метод конечных элементов.
Для дискретизации по времени воспользуемся неявной схемой Эйлера. Пусть АЬ — постоянный шаг по времени. Положим у0 = уо. Для к ^ 1, зная ук-1, находим ук £ И^(О) как решение неравенства
^ У (ук -ук~1)(и-ук)йх + ц У Уук ■ V (и - ук) <1х + \Vu\dx- ! \Чук\йх^ ^ п п п п
^ С(кАЬ) у (и - Ук) йх У и £ Щ.
(6)
Р
Введем обозначения: а = д^, / = С + ау , индекс к в дальнейшем опускаем. Задача (6) эквивалентна (см. [7]) следующей задаче минимизации функционала 3: Ид (О) ^ М: найти функцию у £ Ид (О), такую, что 3(у) ^ 3(и) для любого и £ И^(О), где
</(«) = ^ !(а\и\2 + ц\Чи\2) йх + то J \Vu\dx — ^ /ийх. (7)
п п п
Можно ввести так называемый расширенный лагранжиан, соответствующий задаче (7), следующим образом:
Сг({у, с[}, = ^ У (аМ2 + А1!^!2) йх + то У йх — У /у йх + У А(Уг> — q) йх + ^ J \Уу — q|2 йх, п п п п п
где множитель q имеет смысл градиента скорости, тоА — пластических напряжений.
вестн. моск. ун-та. сер. 1, математика. механика. 2009. № 1
69
Доказано [1], что решение задачи (7) сводится к нахождению седловой точки ({v, q}, Л) G (Hq (Q) х Q) х Q, где Q = (L2(Q))2, функционала Lr, т.е.
Lr({v, q}, X) < Lr({v, q}, Л) < Lr({u, Y}, Л) V{{u, y},X} G (H(Q) х Q) х Q.
С этой целью используем алгоритм типа Узавы из [4] (в литературе он часто называется ALG2).
Для дискретизации вариационных неравенств применяется метод конечных разностей. Обозначим Q = [0, а] х [0,6], h\ = a/Ni и h2 = b/N2. Определим следующие сеточные области: Qi = {xj = (ihi,jh2 ) | i = 0,...,Ni,j = 0,...,N2}, Q2 = {Xij = ((i + 1/2)hi, (j + 1/2)h2 ) | i = 0,...,Ni — 1,j = 0,...,N2 — 1}. Скорость вычисляется в узлах сетки Qi, компоненты тензоров скоростей деформации и напряжений — на сетке Q2. Параметр r, отсутствующий в исходной постановке задачи, является вспомогательным (в алгоритме он фактически выполняет роль итерационного параметра), решение не должно от него зависеть, поэтому для корректности задачи оператор Лапласа аппроксимируется по схеме "косой крест"; для дискретных операторов градиента и векторной дивергенции используется четырехточечный шаблон
vi+i,j+i — vi,j+i + vi+i,j — vi,j vi+i,j+i — vi+i,j + vi,j+i — vi,j
(Vhvh)ij =
2hi 2h2
/ _ А11 + А11 _ А11 \12 _ \ 12 + \ 12 _ \12 \ ' = ^-2Нг-+--) ;
аналогично для (Vh ■ ^н)2 j.
Было бы интересно проверить следующее предположение, которое, к сожалению, доказать пока не удалось, а именно оценить время остановки, используя в (4) вместо 7 другую константу. Кажется привлекательной идея воспользоваться теоретическими результатами для стационарных течений в трубе, которые исследованы намного лучше. В области О могут быть зоны, в которых среда ведет себя как твердое тело (жесткие зоны). Как правило, традиционно рассматриваются два вида жестких зон — зоны застоя, в которых среда покоится, и ядра течения, в которых среда движется как твердое тело. При возрастании то эти зоны увеличиваются, а при достаточно большом то — блокируют течение. В [4] было доказано, что в случае С < 7Тэ\О\-1/2 течение заблокировано (достаточное условие). Заметим, что данное условие совпадает с условием остановки течения для нестационарной задачи (также достаточным). Ранее в [3] было доказано более сильное утверждение: условие С < то/К, где
1 • f п
— = mi
f |Vu| dx
К иеи1(о.)\{о] |Ы|ь1(п) '
является необходимым и достаточным для отсутствия стационарного течения в области О, и приведены значения критического числа для ряда областей (квадрат, прямоугольник, круг, кольцо). Учитывая
М1ь1 (П) < |О|1/2|М|ь2 (П), получим
/\Vu\ йх \О\1/2 / \Vu\ йх
ПП 7 = ш1 ——-- ^ ш!
и€Щ(П)\{о} (П) иеН1(П)\{о> 1и^Ь1 (П)
следовательно,
\ О\1/2
7 < ^ • (8)
Таким образом, в стационарном случае условие С < то7\О\1/2 дает заниженную оценку, а при замене 7 на величину \О\1/2/К получается точная оценка. Предположим, что замена 7 на величину \О\1/2/К в (4), (5) также дает более точные оценки, т.е.
\т\\ь2{п)+< {Шью + «оехр(-^)> «1 = (^г1 -с)(10)
при этом из (8) следует Т\ ^ Тс. Таким образом, представляется возможным явно оценить время остановки. Для квадрата (Ао = 2тг2/а21 К = а/(2 + л/тг), а — сторона квадрата) и для прямоугольника (Ао = 7г2(1 /а2 + 1/Ь2), К = а,Ь/(а. + 6 + (а — Ъ)'2 + ттаЬ), а, Ь — стороны прямоугольника) оценка (9) была проверена численно.
Приведем результаты численных расчетов для следующих данных: при Ь ^ 0 течение стационарное, соответствующее С = 1; при Ь > 0 перепад давления отсутствует, т.е. С = 0. Формы и расположение жестких зон полностью согласуются с теоретическими предсказаниями и оценками [3].
На рисунке приведены графики ||у||ь2 для различных значений то в зависимости от времени (сечение — единичный квадрат). Точки пересечения с горизонтальной осью — время остановки (у = 0). Видно, что жидкость с большим значением то, а следовательно, и с меньшим расходом останавливается быстрее. Пунктирные линии соответствуют результатам численных расчетов, сплошные линии — оценке (10). Разница между вычисленным и оценочным временем остановки составляет менее 5%.
Автор выражает благодарность Л.В. Муравлевой за постановку задачи и полезное обсуждение.
СПИСОК ЛИТЕРАТУРЫ
1. Dean E.J., Glowinski R., Guidoboni G. On the numerical simulation of Bingham viscoplastic flow: old and new results // J. Non-Newtonian Fluid Mech. 2007. 142. 36-62.
2. Георгиевский Д.В. Некоторые неодномерные задачи вязкопластичности: жесткие зоны и устойчивость (обзор) // Изв. РАН. Механ. твердого тела. 2001. № 1. 61-78.
3. Мосолов П.П., Мясников В.П. Вариационные методы в теории жестковязкопластических сред. М.: Изд-во МГУ, 1971.
4. Дюво Г., Лионс Ж.-Л. Неравенства в механике и физике. М.: Наука, 1980.
5. Магомедов О.Б., Победря Б.Е. Некоторые задачи вязкоупругопластического течения // Упругость и неупругость. Вып. 4. М.: Изд-во МГУ, 1975. 152-169.
6. Ильюшин А.А. Деформация вязкопластического тела // Уч. зап. МГУ. Механика. 1940. Вып. 39. 3-81.
7. Гловински Р., Лионс Ж.-Л., Тремольер Р. Численное исследование вариационных неравенств. М.: Мир, 1979.
Поступила в редакцию 09.01.2008
