Численное решение задачи о течении высоковязкой жидкости с монодисперсным наполнителем Текст научной статьи по специальности «Математика»

УДК 519.6:532.529

ЧИСЛЕННОЕ РЕШЕНИЕ ЗАДАЧИ О ТЕЧЕНИИ ВЫСОКОВЯЗКОЙ ЖИДКОСТИ С МОНОДИСПЕРСНЫМ НАПОЛНИТЕЛЕМ

Ю.Н. Константинов, А.А. Ефанов, М.Ю. Альес (Ижевск) Abstract

The numerical research of high-viscous monodispersed suspension flow is presented. The solution is based on Galerkin finite element method for bulk-averaged combined system of equations governing the monophasic continuous media mechanics. 7he solving nonlinear system equations difficulties by means of evolutionary schemes is overcome.

1. Постановка задачи

Рассматривается задача о медленном течении высоковязкой суспензии с низким содержанием монодисперсного наполнителя. Моделирование данного процесса основывается на решении следующей системы уравнений [1]:

div т - grad р - F + p10g = 0, F +(р110 - p!0)al,g = 0, div(u1a1) = 0, div(u*!a1J) = 0,

т = 2ццае, е = 0,5(grad uv + grad ), ца - (11- 1,5а11)/а1, (Г)

uv = и'а1 + u^ot1*,

F = ф(аГ1 )j.i(и1 - и11),

ф(аЛ) = 4,5cT2(l + 1,5(аП )1/3 + 2,25(а1! )2/3 + 1,875(аТ1) + (а11 )4/3) ,

а11 + а1 = 1, pWp10, ри=а7, ап < 0,1,

где и, р, а, j-i, р - соответственно скорость, давление, концентрация, вязкость и плотность; индексы 1 и II соответственно относятся к жидкой и твердой фазе, индексом 0 обозначаются параметры фазы; g - вектор массовых сил; а - диаметр

частиц дисперсной фазы; F- интенсивность межфазного взаимодействия, ф(аП)-функция, учитывающая стесненность обтекания или присутствие соседних частиц.

Система уравнений (1) замыкается следующими граничными и начальными

условиями: и11^ = 0, и11!^ = 0, а11 - а1,1, * = 0.

Дисперсная смесь в данном случае рассматривается как совокупность заполняющих один и тот же объем двух континуумов. Замкнутая система уравнений (1) получена в приближении пространственного осреднения микроуравнений с

помощью ячеечной модели. Вследствие низкой концентрации твердой фазы напряжения в ней не учитываются, а напряжения в жидкой фазе равны средним напряжениям по смеси

Таким образом, решение задачи о медленном течении высоковязкой суспензии

сводится к системе уравнений (1) относительно неизвестных скоростей жидкости и1, твердой фазы и11, давления р и концентраций а!и а!!. Для получения системы уравнений относительно неизвестных скоростей суспензии uv и твердой фазы

и11 выразим из уравнения движения для твердой фазы интенсивность межфазного взаимодействия F и подставим в уравнение движения для жидкой фазы, а также сложим уравнения неразрывности для твердой и жидкой фазы, в результате чего получим следующие две взаимосвязанные системы уравнений:

divx - grad р + g(p! +p[I) = О,

div u v = 0, (2’)

^",\Uv-„V(pIIOV%Vo,

1 - a

div(uI,a1I) = 0. (2”)

Таким образом, моделирование процесса медленного течения при смешении высоковязкой суспензии сводится к решению двух взаимосвязанных систем уравнений для суспензии (2’) и твердой фазы (2”).

2. Метод решения

Рассмотрим для сравнения методы решения систем уравнений (1) и (2) Обе системы уравнений решаются методом конечных элементов в слабой формулировке Галеркина [2]. Аппроксимация производится с помощью треугольных конечных элементов. После применения метода взвешенных невязок система уравнений (1) будет выглядеть следующим образом:

| т grad WdV + Jр div WdV + j FWafK - f pIOgWdV = 0,

a a n a

Jf\VJK +jw(pIIO-pl°)aI1grfK = 0, (3)

q a

{¥ divlVa1 )dV = 0, JVdivfiiVyr = 0,

П Q

где W, *P - кусочно-линейные весовые функции.

Для данных задач распределения скоростей, давлений и концентраций аппроксимируются комбинацией базисных функций Ne и степеней свободы решения, в качестве которых для лагранжевых конечных элементов выступают узловые значения му ,е , и]е, ul}e ,ale,aI1 ,ре неизвестных функций:

а1 = а\Nе, аГІ = а|,1Л,е, р = реИе,

(4)

Ме ~ У'1 + У/ ех1 . е~\^, / = 1,2 ; т- 1, А,

где V,, у,е - коэффициенты, зависящие от координат узлов [2], Ь - общее число узлов.

Применение процедуры метода конечных элементов с использованием слабой формы Галеркина и аппроксимаций (4) к системам уравнений (1) и (2) приводит для

(1) к системе нелинейных алгебраических уравнений (относительно компонент

скорости жидкой и1 и твердой фазы и11, концентраций а1 и а11, давления р, для (2) к двум взаимосвязанным системам проекционно-сеточных уравнений (относительно компонент скорости суспензии иу и твердой фазы и11, концентрации а11, давления р). Во втором случае система уравнений для твердой фазы является нелинейной.

Для решения системы нелинейных алгебраических уравнений воспользуемся методом [3], который основан на расширении исходной граничной задачи и сведении ее к задаче Коши. Для этого в системах уравнений (1) и (2”) уравнения неразрывности для твердой фазы заменяются эволюционным соотношением:

г де £, - искусственное время Используя явную схему с шагом Ас для аппроксимации (5), получим следующую итерационную модель:

где / - номер шага по А£,. Анализ схемы (6) показал, что сходимость достигается при следующем условии:

Поскольку для решения системы нелинейных алгебраических уравнений используется итерационный метод, то для преодоления сложностей, связанных с расчетом давления и удовлетворением условию неразрывности, предпочтительнее использовать метод искусственной сжимаемости [4], в котором уравнение неразрывности суспензии заменяется эволюционным соотношением:

(5)

а11 - а11-А^Іу(иПаП),

(6)

А^ > О, с!іу(и11 О> 0 > где £ - ошибка на итерации.

Для реализации условия (^(и11 0> 0 выражение (6) следует решать в виде. (+1 / | /

(6;)

(7)

Уравнение (7) может быть аппроксимировано по неявной схеме с шагом Аб,

(+1 I

аі;

Из соотношения (8) может быть выражено давление

(8)

/+1 ; 1+\

р -- р- ХсИу и (9)

и подставлено в уравнение движения с целью понижения числа уравнений системы проекционно-сеточных уравнений. Параметр X необходимо выбирать таким образом, чтобы обеспечить сходимость к стационарному численному решению уравнений Стокса. Сходимость достигается при следующем условии [5]:

.24

и 2

h V

А^

Re

< 1, СҐ > О,

(10)

где 11е - число Рейнольдса; И - характерный шаг сетки.

Системы алгебраических уравнений, полученные в результате применения к (1) и

(2) процедуры метода конечных элементов с использованием слабой формы Галеркина, аппроксимаций (4), а также итерационных схем (6) и (9), в операторном виде будут выглядеть следующим образом.

(

А'(аІІ) + М(аІІ) + Да11)

f+1

«Л

Л/(аИ) +А/(аИ) + Да11)

г + 1 II

и11 = К(р) + G ,

;+1

( + 1

А(ап)и1 + Л(а11)и1Г = ^'(аП).

01)

В а

/+1

.1

; + 1 11

/ і +1

а* = 1 - а н-1 I Р

1+1

11

j I +1 /' +1 / +1 І +1

S(p)~ D(ul ,ип, а1, а11)

if 1 %

Nuv =G(all)-K(p), (12’)

A p - S(p)~ D{uv), g j g

L uu = M(ull,uv,ctU), 0Г')

я+i g g

Ba{1 = С(аи,ип).

/ i Здесь N(ali) - матрица, содержащая члены с вязкостью; М(а11) - матрица,

/

включающая параметр искусственной сжимаемости, /-(а11) - матрица, содержащая

/

интенсивность межфазных взаимодействий в суспензии; К{р)- вектор, включающий в

!

себя давление с предыдущей итерации; F(aJ1) - вектор массовых сил.

Таким образом, исходная задача на классе сеточных функций сведена в первом случае к системе уравнений (11) относительно компонент скорости жидкой и1 и

твердой фазы и11, концентрации а1 ,а^и давления р. Система уравнений движения

решается при фиксированном распределении концентрации а11, система уравнений неразрывности твердой фазы решается при фиксированном распределении скорости

и11, система уравнений для давления решается при фиксированном распределении скоростей и1, и11 и концентраций а11, а1 Во втором случае задача сведена к двум взаимосвязанным системам уравнений (12), (12”) относительно скоростей иу,и^, давления р и концентрации а11. Система уравнений для суспензии и твердой фазы решается при фиксированных распределениях концентрации а11 и скорости иу соответственно

Полученные из уравнений движения в первом (11) и во втором случае (12) системы алгебраических уравнений имеют различные матрицы коэффициентов В первом случае это несимметричная матрица, а во втором - симметричная. Это обстоятельство предполагает использование различных методов решения систем алгебраических уравнений. В первом случае применяется метод Гаусса, во втором метод Холесского. Известно, что метод Холесского в отношении потребляемых ресурсов компьютера намного выгоднее метода Гаусса, в частности, для решения системы линейных алгебраических уравнений с т неизвестными при использовании

метода Гаусса потребовалось бы - гг? вычислительных операций, в то время как при

использовании метода Холесского - 2/и2 (6] В первом случае приходится хранить в памяти ЭВМ всю матрицу и перестраивать ее на каждой итерации, что невыгодно по сравнению со вторым алгоритмом (12), в котором матрицы на временном шаге создаются только один раз, а перестраивать на каждой итерации приходится только правые части систем уравнений. Все преимущества второго алгоритма (12) позволяют значительно ускорить получение решения системы уравнений, описывающей течение суспензии, а также дают возможность получать решение при использовании

одинакового с первым алгоритмом размера памяти, но на более мелкой сетке. В частности, для выполнения расчетов одного временного шага течения суспензии на сетке с числом узлов 121 на компьютере с процессором 1486 потребовалось: с помощью первого алгоритма ~ 10 секунд, а с помощью второго алгоритма - 1 секунда

Переход на новый временной шаг осуществляется при использовании

1+1 I

следующего кинематического соотношения: х; = х14 и| Дг, где Д1- шаг по времени, который выбирается из условия Д1«Ь/2,Ь - шаг сетки.

Проведенные численные исследования [7] при использовании алгоритмов (11) и (12) показали полную идентичность полученных результатов, а также подтвердили значительное превосходство в эффективности второго алгоритма над первым

Таким образом, рассмотрены два алгоритма решения системы уравнений,

описывающей медленное течение высоковязкой суспензии Показана возможность эффективного применения эволюционных методов для решения данной задачи,

которые имеют значительные преимущества по сравнению с другими методами.

Работа выполнена при поддержке РФФИ (код проекта 95-01-00293а) и 158ЕР (грант № 520с1)*

Библиографический список

1. Нигматуллин Р.И. Основы механики гетерогенных сред - М .: Наука, 1973 -336с.

2. Коннор Д., Бреббиа К. Метод конечных элементов в механике жидкости. - Л : Судостроение, 1979. - 263 с,

3. Липанов А.М., Альес М.Ю, Константинов ЮН. Численное моделирование ползущих течений неньютоновских жидкостей со свободной поверхностью// Математическое моделирование, - 1993. - Т.5 - №7. - С.3-9

4. Темам Р. Уравнения Навье-Стокса Теория и численный анализ - М.: Мир, 1981. - 408 с.

5. Пейре Р., Тейлор Т. Вычислительные методы в задачах механики жидкости - Л.: Гидрометеоиздат, 1986 - 352 с.

6 Амосов А.А., Дубинский Ю.А., Копченова Н.В. Вычислительные методы для инженеров. Учеб. пособие. - М : Высш шк., 1994. - 544. с.

7. Константинов ЮН, Ефанов А А., Альес М.Ю. Численные исследования смешения монодисперсной суспензии// Сб. трудов конференции “Применение математического моделирования для решения задач в науке и технике”, -Ижевск, 1996.-С.227-233.