CC BY
13тогда для определения p(x) получим квадратное уравнение (а*/(1 — ß2))z2 + Coz — x = 0 z = pßl] параметры а* и Co находятся из граничных условий. Решение записывается в форме
V = f1(x,x0,t)=zW (а+201^2)-
Во втором случае положим wi = —а*ßipßl/2+1 /(ß\ + 2), аналогично предыдущему получим C0z2 - a*z/( 1 - ß2) - х = 0, г = р^2, у = h(x, x0,t) = z2'^ (^i2x - .
Анализ этих решений не проводится; здесь они указаны как возможные варианты дальнейших исследований.
Замечание 3. В предложенной работе случаи анизотропии (фактически ортотропии) свойств материала и контактного трения рассмотрены отдельно; обобщение на совместное действие обоих эффектов при условии совпадения главных осей ортотропии дано в публикации [5], где показано, что по постановке задача не отличается от рассмотренных. Задача о течении в слое при несовпадении главных осей ортотропии требует специального исследования.
СПИСОК ЛИТЕРАТУРЫ
1. Ильюшин A.A. Вопросы теории течения пластического вещества по поверхностям // Прикл. матем. и механ. 1954. 18, вып. 4. 265-288.
2. Ильюшин A.A. Полная пластичность в процессах течения между жесткими поверхностями, аналогия с песчаной насыпью и некоторые приложения // Прикл. матем. и механ. 1955. 19, вып. 6. 693-713.
3. Кийко И.А. Пластическое течение металлов // Научные основы прогрессивной техники и технологии / Под ред. В.П. Лымзина. М.: Машиностроение, 1985. 102-133.
4. Кийко И.А. О форме пластического слоя, сжимаемого параллельными плоскостями // Прикл. матем. и механ. 2011. 75, вып. 1. 14-25.
5. Кийко И.А. Анизотропия в процессах течения тонкого пластического слоя // Прикл. матем. и механ. 2006. 70, вып. 2. 344-351.
6. Кийко И.А. Вариационный принцип в задачах течения тонкого слоя пластического материала // Докл. АН СССР. 1964. 157, № 3. 551-553.
7. Полянин А.Д., Зайцев В.Ф. Справочник по нелинейным уравнениям математической физики. М.: Физматгиз, 2002.
8. Безухое В.И. Об осадке пластического слоя некруговой формы в плане: Канд. дне. М., 1955.
9. Зайцев В.Ф., Полянин А.Д. Справочник по обыкновенным дифференциальным уравнениям. М.: Физматлит, 2001.
10. Степанов В.В. Курс дифференциальных уравнений. М.: Физматгиз, 1959.
Поступила в редакцию 24.06.2013
УДК 532.685
ТОЧНЫЕ РЕШЕНИЯ НЕЛИНЕЙНЫХ УРАВНЕНИЙ ТЕЧЕНИЯ СУСПЕНЗИИ В ПОРИСТОЙ СРЕДЕ
Н. Е. Леонтьев1, Д. А. Татаренкова2
Построено общее решение нелинейной системы уравнений, описывающей одномерные течения малоконцентрированной суспензии в пористой среде с учетом оседания частиц. Указан ряд решений, выражающихся в элементарных функциях. Проанализированы условия образования нефизичных особенностей.
1 Леонтьев Николай Евгеньевич — канд. физ.-мат. наук, доцент каф. гидромеханики мех.-мат. ф-та МГУ, e-mail: leontiev_nQmail.ru.
2 Татаренкова Дарья Александровна — студ. каф. гидромеханики мех.-мат. ф-та МГУ, e-mail: [email protected].
Ключевые слова: течение в пористой среде, суспензия, фильтрация, отложение частиц, градиентная катастрофа.
The general solution to a nonlinear system of differential equations describing one-dimensional flows of a low-concentrated suspension through a porous medium with consideration of particle deposition is constructed. A number of closed-form solutions are specified in terms of elementary functions. The formation of non-physical singularities is studied.
Key words: flow through a porous medium, suspension, deep bed filtration, particle deposition, gradient catastrophe.
1. Рассматривается течение малоконцентрированной суспензии в пористой среде с учетом оседания частиц на пористый скелет. При обычных предположениях [1, 2] объемная концентрация взвешенных частиц а, пористость m, давление p и скорость фильтрации суспензии u находятся из системы, содержащей уравнения баланса массы для взвешенных частиц и всей суспензии, уравнение баланса импульса (закон Дарси), а также наиболее употребительный [3] вариант кинетического уравнения
д(та) , , . dm , п(а) dm .. . . .
— V div(cm) = divu = 0, 0 = - gradp - щyu, — = -7/(m)a|u|,
где t — время; y — константа, обратная величина к которой задает характерный пространственный масштаб задачи; f (т) — известная положительная монотонная функция, вид которой зависит от механизма удержания частиц на скелете; п(а) — динамическая вязкость суспензии; k(m) — проницаемость пористой среды.
Для одномерных течений вдоль оси x с плоскими волнами без ограничения общности [2] можно считать модуль скорости фильтрации uo постоянным. В этом случае после введения безразмерных времени t1 = uoYt и координаты X = yx получается нелинейная система (далее штрихи над безразмерными переменными опускаются)
д(та) да dm дт „. . . .
уравнение баланса импульса решается отдельно и служит только для нахождения распределения давления в суспензии. В практических приложениях для исследования системы (1) (или аналогичных ей в случае учета других факторов) обычно применяются численные методы или используются двусторонние оценки [3, 4], однако известны и относительно немногочисленные аналитические решения подобных задач (укажем решения начально-краевой задачи [1], автомодельные решения [5] и решения в виде бегущей
(1)
2. Исключив а с помощью второго уравнения в (1), получаем уравнение для пористости m(x,t):
di di 1 т {т9(т)^) + д:х {9{m)^J +Ж = 0> д{т) = т>0, (2)
которое, как можно заметить, однократно интегрируется:
dm> dm
mg(m) ——Ь g(jn) —--Ь т + Ф(ж) = 0. (3)
dt dx
В последующем мы рассмотрим задачу Коши для системы (1) когда в момент t = 0 на всей числовой оси задаются начальные пористость mo(x) и концентрация ао(х), что определяет постоянную интегрирования Ф:
Ф(ж) = т0а0 - gimo) ^^ - ш0.
Такая постановка задачи — течение в бесконечном пласте — является модельной, так как, во-первых, в этом случае для осуществления течения формально требуется бесконечный перепад давления (аналогично, например, течению Пуазейля в бесконечной трубе), а во-вторых, в практических приложениях пористый образец всегда конечен. С другой стороны, это ограничение несколько смягчается тем обстоятельством, что при сужении решения на конечный интервал будет получаться некоторое решение начально-краевой задачи для конечного пласта, соответствующее специальному выбору концентрации на одной из его границ.
Интегрирование квазилинейного уравнения (3) стандартным образом сводится к решению автономной системы обыкновенных дифференциальных уравнений (ОДУ)
Щ^=д{т{,з)), <^- = т(8)д(т(8)), = -Ф(ф)) - ш(в);
ж(0) = Хо, т=0, т(0)= то (хо),
где в — параметр вдоль характеристик, хо — пространственная координата точки на характеристике в момент в = 0 (идентификатор конкретной характеристики). Отметим, что начальное распределение концентрации ао учитывается через функцию Ф(ж).
Дальнейшие рассуждения зависят от вида функции д(т) (в терминах исходной системы уравнений — от вида функции f (т)).
Если д(т) строго монотонна, так что существует обратная функция д-\ то для функции Е(ж), входящей в параметрические зависимости
йж(в) сРж(в) йЕ (ж)
получим задачу Коши для нелинейного ОДУ первого порядка:
В этом случае общее решение уравнения (2) представляется в параметрическом виде
х
ж = ж(х, Хо) = Х, I = Кх,Хо ) = ! С(хг, Хо) Х т = т(х,Хо) = С(х,Хо), С(х,Хо) = д-1(Е (х,Хо)),
хо
вХ
в обозначении функции Е явно подчеркнута з^исимость от параметра Хо, указанного в качестве второго аргумента. Наконец, для концентрации а с учетом тождества = 0 получим выражение
дт дС , .
(X) Хо)
а = а(х,х о) = "% =-%
91 то(хо) — [ тг—(хъ Хо)
дХо хо дХо
В специальном случае постоянной функции д(т) (можно ститать, что д(т) = 1), который важен с практической точки зрения, обратная функция не существует, поэтому решение системы (4) записывается в другом параметрическом представлении:
(х \ х
то(хо)еХ0 - ! Ф(Х1 )ех1 йхА , í = ¿(х,Хо) = У т(х1,Хо) йхг,
хо хо
а распределение концентрации определяется выражением
ао(хо)ехр(хо - х)
а = а(Х, Хо) =
1 + ао (хо)(ехр(хо - х) - 1)'
Из последней формулы вытекает явный критерий наступления градиентной катастрофы (нефизичного обращения концентрации а в бесконечность): если на всей числовой оси начальная концентрация ао (ж) < 1, то во всех конечных точках концентрация а(ж,Ь) не обращается в бесконечность; если же в начальный момент имеется точка, в которой ао > 1, то в какой-то момент концентрация обязательно обращается в бесконечность.
Во избежание недоразумений специально подчеркнем, что, во-первых, здесь мы не учитываем физических условий малости концентрации (а ^ 1) и конечности пористости (т < 1), выполнение которых, естественно, должно дополнительно проверяться, а во-вторых, система уравнений в случае концентрации,
формально превосходящей единицу (или концентрацию, соответствующую наиболее плотной упаковке взвешенных частиц) в какой-либо области, может описывать реальные процессы, если мы интересуемся течением в другой пространственной области, где все естественные физические ограничения выполняются (так, например, происходит в приводимом ниже примере). Отметим также, что при резком росте концентрации рассматриваемые модели перестают быть физически применимыми и должны заменяться на более полные [1].
Если засорение пористого пространства незначительно (m(x,t) ~ mo(x)), то вместо системы (1) часто используют ее упрощенный вариант [2]:
, . да да dm dm .. .
+ ^ = -m=-fim)a> (6)
в котором кинетическое уравнение не подвергается упрощению (допускается значительное изменение функции f (m)). В этом случае решение начальной задачи может быть найдено аналогичным образом:
х
х = х(х,Хо) = X, m = т(х,Хо), t = t(x,Xo) = j mo(xi) dxi,
X0
где m(x, Xo) находится из решения задачи Коши для ОДУ
dm(x) Ф(х) + m(x)
dx g(m(x))
m(xo) = mo(xo), X > Xo•
3. Приведенные выше формулы для неупрощенной системы (1) в ряде случаев позволяют найти выражающиеся в элементарных функциях решения, которые удобно использовать, например, при тестировании численных методов.
Для случая д(т) = 1 отметим прежде всего единственное автомодельное решение т = —Ь, а =
1
переменных (с1, С2 — константы):
С2 — Ь 1 т = ----, а =
1 + в\ exp(x)' 1 + в\ exp(x)
при m0(x) = c2/(1 + c\ exp(x)), a0(x) = 1/(1 + c\ exp(x));
C2 exp(—x)
m = exp(—x)\/ci — c2t, a =
Vci - c2t
при шо(ж) = у/сГехр(—х), ао(х) = ехр(—х)/{2^/с\).
Укажем примеры более громоздких решений в случае g(m) = 1:
Ci(l - с2) In (l + 171 ~СЛ +1 + rn - а = 0, ci(c2 - 1)(1 - а) In ^ ~ + (t - а)а - t + cic2 = 0 V CiС2 J C2(a — 1)
при m0(x) = ci, a0(x) = c2;
(ci(l — c2) — t — m) ехр(ж) + (C\C2 —t — m) exp | шехР(ж)—ж j = q,
V ciC2 J
—c\c2a In (——1 — (ci + tехр(ж))а + C\C2 = 0 \c2(a - 1)/
при mo(x) = c1 exp(—x), ao(x) = c2;
C1(C2 - 1) In (exp + ln(1 + ехр(ж))) +
+ (c1 — m) exp(x) + c1 (1 — c2 )x — t — m + c1 =0,
С2(аexp(x) + а — 1) а(с2 - 1)
-cic2(exp(x) +1) In а (т) + (а_1) + ci(c2 " 1) ехр(.т) In _ +
- 2 - 2
+ cic2(1 + exp(x)) ln(1 + exp(x)) — t exp(x) = 0
при m0(x) = Сь a0(x) = c2/(1 + exp(x));
, . , . 2 / (c1 — m) exp(x) \ 2 , .
(ci — m) ехр(ж) — c2t — c\c2 exp--x + c\c2 exp(— x) — с2т + cic2 = 0,
V c1c2 )
ac1 (1 + c2 exp(—x)) ln при m0(x) = ci: a0(x) = c2 exp(—x).
c2a
a exp(x) + c2(a — 1)
+ (t + c1)a — c1c2 exp(—x) = 0
(7)
(8)
Рис. 1. Пористость (7) и концентрация (8) для течения в конечном пласте при е\ = 0,1, С2 = 0,1; кривая на плоскости (ж,£) — точки с бесконечными производными
Для последнего решения на рис. 1 показан вид функций т(х, ^ и а(х, ¿) на конечном участке пласта при С1 = 0,1 е2 = 0,^. В этом примере начальная концентрация ао на рассматриваемом интервале не превосходит единицы, однако на всей числовой оси ао может превышать единицу, что и объясняет наличие градиентной катастрофы. Отметим качественное отличие приводенного решения от решения для упрощенной модели (6)
m(x, t) = c1 — c2te
a(x, t) = c2e
при тех же начальных условиях, для которого образование бесконечных производных не происходит никогда. Это различие связано с тем, что для (6) характеристики на плоскости (x,t), определяемые до решения задачи, не имеют пересечений.
При непостоянной функции g(m) в ряде случаев уравнение (5) также может быть проинтегрировано в конечном виде. Например, если начальные данные подобраны так, что $(x) = const, то оно интегрируется в квадратурах. В качестве иллюстрации на рис. 2 приводится вид
m(x, t) (1) (6)
лей в случае g(m) = k1 ln m + k2 при начальных данных m0(x) = c1 exp(—x), a0(x) = k1(x — lnc1) — k2 + 1, для которых $(x) = 0 (в этом случае решение выражается в элементарных функциях, но очень громоздко).
4. При g(m) = 1 имеется аналогия между одномерными течениями с плоскими и цилиндрическими волнами, что позволяет переносить полученные выше решения на оеееимметричеекий случай. В самом деле, оеееимметричеекие течения описываются системой
Рис. 2. Сравнение решений для полной моде-(1)
ли (6) (белая заливка); е\ = 0,1 к\ = —0,1, к2 = 0,5
д (ma) 1 d(rua) dm dm dt r dr dt ' dt
= —^а\и\, u(r, t) =
q(t)
r
где r — расстояние до оси симметрии; q(t) — функция, пропорциональная объемному расходу суспензии (в расчете на единицу длины вдоль оси симметрии); u — радиальная составляющая скорости фильтрации. После введения новой переменной a(r,t) = r • m(r, t) система приобретает вид
д (да) да да да
= w = q = q(i)'
формально совпадающий с видом системы для течений с плоскими волнами, однако теперь простран-r
течений при g(m) = 1 для тех точек на плоскости (r,t), область зависимости [7] которых содержится в интервале r > 0, решение очевидным образом переносится со случая течений с плоскими волнами.
Работа выполнена при поддержке РФФИ (грант № 14-01-00056).
СПИСОК ЛИТЕРАТУРЫ
1. Леонтьев Н.Е. О структуре фронта пористости при движении суспензии в пористой среде // Вестн. Моск. ун-та. Матем. Механ. 2006. № 5. 73-76.
2. Леонтьев Н.Е. Засорение пористого пласта при движении фронтов с конечным скачком пористости // Вестн. Моск. ун-та. Матем. Механ. 2009. № 5. 75-80.
3. Dqbrowski W. Consequences of the mass balance simplification in modelling deep filtration // Water Res. 1988. 22, N 10. 1219-1227.
4. Dqbrowski W. A method of calculating the concentration of deposit around an injection well // Water Res. 1984. 18, N 6. 709-717.
5. Кузнецов А.Ю., Поелавекий С.А. Исследование математической модели механической суффозии // Вестн. Харьк. нац. ун-та. Матем., прикл. матем. и механ. 2009. № 875. 57-68.
6. Кон С., Леддер Г., Логан Д. Анализ модели фильтрации в пористой среде // Матем. модел. 2001. 13, № 2. 110-116.
7. Курант, Р. Уравнения с частными производными. М.: Мир, 1964.
Поступила в редакцию 14.02.2014
